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The exterior gravitational field of a slowly-rotating neutron star can be characterized by its multi- 
pole moments, the first few being the neutron star mass, moment of inertia, and quadrupole moment 
to quadratic order in spin. In principle, all of these quantities depend on the neutron star's internal 
structure, and thus, on unknown nuclear physics at supra-nuclear energy densities, all of which is 
usually parameterized through an equation of state. We here find relations between the moment 
of inertia, the Love numbers and the quadrupole moment (I-Love-Q relations) that do not depend 
sensitively on the neutron star's internal structure. Such universality may arise for two reasons: (i) 
these relations depend most sensitively on the internal structure far from the core, where all realistic 
equations of state mostly approach each other; (ii) as the NS compactness increases, the I-Love-Q 
trio approaches that of a BH, which does not have an internal-structure dependence. Three im- 
portant consequences derive from these I-Love-Q relations. On an observational astrophysics front, 
the measurement of a single member of the I-Love-Q trio would automatically provide information 
about the other two, even when the latter may not be observationally accessible. On a gravitational 
wave front, the I-Love-Q relations break the degeneracy between the quadrupole moment and the 
neutron-star spins in binary inspiral waveforms, allowing second-generation ground-based detectors 
to determine the (dimensionless) averaged spin to 0.01, given a sufficiently large signal-to-noise ratio 
detection. On a fundamental physics front, the I-Love-Q relations allow for tests of General Rela- 
tivity in the neutron-star strong-field that are both theory- and internal structure-independent. As 
an example, by combining gravitational-wave and electromagnetic observations, one may constrain 
dynamical Chern-Simons gravity in the future by more than 6 orders of magnitude more stringently 
than Solar System and table-top constraints. 

PACS numbers: 04.30.Db,97.60. Jd 
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I. INTRODUCTION 

Neutron-star (NS) astrophysics can provide crucial in- 
formation about nuclear, gravitational-wave (GW) and 
fundamental physics that would be difficult to obtain by 
other means. On a nuclear physics front, NS observations 
allow us to probe the equation of state (EoS) of nuclear 
matter [1] well beyond the densities available in Earth 
laboratories. For example, observations of the NS mass- 
radius relation and the mass-moment-of-inertia relation 
can be used to infer NS EoS within a certain observa- 
tional uncertainty [2-7]. 

On a GW physics front, the detection of GWs emitted 
during the late inspiral and merger of NS binaries could 
also be used to extract information about the EoS [8- 
10]. Binary NSs are, in fact, one of the most promising 
GW sources [11-13] for second-generation, ground-based 
detectors, such as Adv. LIGO [14], Adv. VIRGO [15] 
and KAGRA [16]. Since NSs are tidally deformed in 
the late inspiral and merger, violating the test-particle 
approximation, NS binary waveforms must include cor- 
rections induced by the NS internal structure, for exam- 
ple in terms of NS tidal Love number [8-10]. In turn, 
this implies that a sufficiently high signal-to-noise ratio 
(SNR) detection of such a GW could be used to extract 
information about the NS EoS [17-30]. 

On a fundamental physics front, NSs are ideal to test 
General Relativity (GR) since they produce strong gravi- 
tational fields [31-33]. Currently, GR has passed all Solar 
System tests with flying colors, but these only sample the 
weak field regime [34, 35], where gravitational fields are 
stationary and weak, and all characteristic velocities are 
much smaller than the speed of light. Electromagnetic bi- 
nary pulsar observations can test GR in a regime where 
the gravitational field is much stronger than in the Solar 
System, but still sufficiently non-dynamical that one can 
expand in the ratio of the orbital velocity to the speed of 
light to leading-order [31]. 

The exterior gravitational field of NSs, however, is not 
just determined by their mass and radius, but also by 
higher multipole moments, like the moment of inertia and 
the quadrupole moment, and ignorance of the NS EoS 
can hinder the extraction of the latter from observations. 
On the GW physics front, degeneracies between the NS 
spin and the quadrupole moment prevent future detec- 
tions from separately measuring these quantities. On the 
fundamental physics front, degeneracies between the ef- 
fect of the NS EoS and modified gravity corrections on 
observables prevent robust tests of GR that are internal- 
structure independent. 

In this paper, we take the first steps toward resolving 
this problem by discovering new relations between the 
NS moment of inertia /, the NS Love numbers and the 
(spin-induced) NS quadrupole moment Q^''"*) (I-Love-Q 
relations) that are essentially EoS independent for slowly- 
rotating NSs [36]. Physically, the moment of inertia 
quantifies how fast a NS can spin given a fixed spin angu- 
lar momentum 5, the quadrupole moment describes how 



much a NS is deformed away from sphericity, and the 
Love number characterizes how easy or difficult it would 
be to deform a NS. 

The moment of inertia. Love numbers and quadrupole 
moment can be computed by numerically solving for the 
interior and exterior gravitational field of a NS in a slow- 
rotation [37] and a small tidal deformation approxima- 
tion [8] , to quadratic order in the former and to linear or- 
der in the latter. The moment of inertia and quadrupole 
moment can be obtained from the asymptotic behavior 
of the (t, (f)) and (t, t) components of the metric at spa- 
tial infinity respectively, which depend on the interior 
solution through matching boundary conditions at the 
NS surface that ensure metric continuity and differentia- 
bility. Although the moment of inertia is a first-order 
in spin quantity, the quadrupole moment is generated by 
quadratic spin terms. The tidal Love number A^*"^-* is de- 
fined by the ratio between the tidally-induced quadrupole 
moment and the tidal field due to a companion NS, which 
can be calculated in a similar fashion. 

One would expect that all of these quantities should de- 
pend quite sensitively on the NS EoS; after all, a fluffier 
star should be more easily deformablc than a stiffer star. 
We flnd here, however, that, these quantities seem to sat- 
isfy almost universal relations when plotted against each 
other that are essentially independent of the NS EoS. 
Figure 1 shows the I-Love (left) and Q-Love (right) re- 
lations, where /, A^'"^^ and Q^™') are normalized to M^, 
Ml and Ml{S/M'^f respectively, with the NS mass 
and S its spin angular momentum. The different curves 
represent the relations using different EoSs (APR [38], 
SLy [39], Lattimer-Swesty (LS220) [40], Shen [41, 42] 
and polytropic EoSs with indices of n = 0.6, 0.8 and 
1). The symbols represent numerical solutions, while the 
solid curve is a single fitting function. The bottom of 
this figure shows the fractional errors between the fitting 
function and the numerical results. Observe that these 
relations are EoS independent to within 0(1)%. 

We have found two possible reasons that could explain 
such a weak EoS dependence. The first is that the I-Love- 
Q trio may depend most sensitively on the NS outer lay- 
ers, far from the core, where all realistic EoSs approach 
each other. In this interpretation, the TLove-Q relations 
do depend on the EoS, but only in a regime where the 
EoSs contributes similarly to the I-Love-Q trio. The sec- 
ond reason is based on the fact that the I-Love-Q trio for 
NSs approaches the TLove-Q relations for a BH, as one 
increase the NS compactness. For BHs, these relations 
are clearly independent of the BH internal-structure (or 
lack thereof) due to the the no-hair theorems [43-48], 
which lead to well-known expressions for all multipole 
moments in terms of just the mass and spin [49, 50]. But 
for NSs, such expressions do not exist because there is no 
NS no-hair theorem. In spite of this, we still find a NS 
universal relation between the moment of inertia (and 
thus the spin angular momentum) and the quadrupole 
moment, similar to that which arise for BHs due to the 
no-hair theorems. 



3 




FIG. 1. (Color Online) (Top) Fitting curves (solid curve), given in Eq. (54), and numerical results (points) of the universal I-Love 
(left) and Q-Love (right) relations for various EoSs. These quantities are normalized as follows / = I/M't, A'"'*' = X'^'-'"^^ / Ml 
and Q = Q'^'"'^'' /\Ml{S/Mlf]. The parameter varied along each curve is the NS central density, or equivalently the NS 
compactness, with the latter increasing to the left of the plots. For reference, we also show the corresponding NS mass for the 
APR EoS on the top axes and a vertical dashed line when Af, = IMq. (Bottom) Relative fractional errors between the fitting 
curve and numerical results. Observe that these relations are essentially EoS independent, with loss of universality at the 1% 
level. 



The universal I-Lovc-Q relations tell us that there is 
an effacing of internal structure in play here, i.e. the 
expected internal-structure dependence of the I-Love-Q 
relations effaces away. One might think that such an 
effect is a consequence of the celebrated effacement prin- 
ciple [51] in GR. Hovifever, this is not quite right be- 
cause the latter states that the motion of compact ob- 
jects is independent of their internal structure; the ef- 
facement principle says nothing about the multipolar- 
decomposition of the object's gravitational field or of its 
tidal deformations. Of course, the effacement principle 
holds in GR for BHs, but it is violated for NSs, with 
internal-structure corrections to the center of mass ac- 
celeration entering first at 5 post-Newtonian (PN) or- 
der^ for systems of non-spinning bodies. On the other 
hand, the I-Love-Q relations interconnect different mul- 
tipolc components of the exterior gravitational field of 
isolated bodies, saying nothing about their relative mo- 
tion. 

The I-Love-Q relations have immediate applications 
to observational astrophysics, GWs and fundamental 
physics, breaking degeneracies that would otherwise pre- 
vent us from taking full advantage of NS observations. 
On the observational astrophysics front, the measure- 



^ A term of Ath PN order is suppressed relative to the leading- 
order term by a factor of 0{v'^'^ / c^^), where v is the character- 
istic velocity of the system and c is the speed of light. 



ment of any single member of the I-Love-Q trio would au- 
tomatically provide information about the other 2 mem- 
bers, even if the latter arc not easily accessible from an 
observational viewpoint. For example, if one could mea- 
sure the moment of inertia of the primary NS of the dou- 
ble binary pulsar J0737-3039 [52-54], one could then ob- 
tain its quadrupole moment and its tidal Love number 
through the I-Love-Q relations without any further mea- 
surements. This is particularly important because the 
Love number and the quadrupole moment cannot be eas- 
ily extracted from binary pulsar observations, since they 
have a much weaker effect on observables. 

On the GW physics front, the I-Love-Q relations can 
break the degeneracy between the NS quadrupole mo- 
ment and the NS spins, given a sufficiently large SNR 
detection of a NS binary inspiral. The first spin-induced 
modification to the waveform, a spin-orbit coupling, en- 
ters at 1.5 PN order in the waveform phase [55]. Given 
a large SNR detection, one can then extract this phase 
term, and thus measure a certain combination of the in- 
dividual spins. In order to extract both spins, however, 
one needs to also measure the spin-spin correction to the 
waveform, which enters at 2PN order. At this same or- 
der, however, the quadrupole moment also modifies the 
waveform phase, leading to a 100% degeneracy between 
Qlrot) and the individual spins. 

The Q-Love relation can be used to break this degener- 
acy. One can write the quadrupole moment as a function 
of the Love number, which enters at 5PN order in the 
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FIG. 2. (Color Online) Measurement accuracy of spin pa- 
rameters /3, Xa and Xa with Adv. LIGO given a detection at 
a luminosity distance of 100 Mpc with SNR ~ 30. We con- 
sider three different NS binaries, labeled by (i), (ii) and (iii), 
as described in the text. A/? and Axa are almost identical 
among all systems, since they are dominated by their priors 
(|/3| < 0.2 and \xa\ < 0.1). Observe that Xa can be measured 
to approximately 0.01 when we use the Q-Love relation to 
break spin degeneracies. 



waveform phase [17]. This forces a correlation between 
the quadrupole moment piece of the 2PN term and a 
5PN term that is weakly correlated with other binary 
parameters. Recently, [25, 26] suggested that second- 
generation, ground-based detectors could be used to ex- 
tract the Love number. Therefore, such measurement of 
the tidal Love number, in combination with the Q-Lovc 
relation, determines the NS quadrupole moment, which 
then allows for a measurement of the averaged spin pa- 
rameter Xs = (xi + X2)/2, where xa is the individual 
(dimensionless) spin parameter of NS A. 

Figure 2 shows the projected measurement accura- 
cies of spin parameters Xa = (xi ~ X"!)!^! Xs and the 
1.5PN phase term /? as functions of xi for 3 differ- 
ent systems: (i) (mi, 7712) = (1.45, 1.35)Mo, xi = X2, 
(ii) (TOi,m2) = (1.45,1.35)Mo, Xi = 2x2 and (iii) 
(mi, 771,2) = (1.4, 1.35)Mq, Xi = X2, where wia is 
the NS mass of the A component. We used second- 
generation, ground-based detectors and luminosity dis- 
tance of lOOMpc for SNRs ^ 30. The measurement ac- 
curacies of /3 and Xa are almost identical among the 3 
systems because they arc dominated by their priors. On 
the other hand, one can measure the averaged spin Xs 
to 0.01 if one uses the Q-Lovc relation. Such a measure- 
ment accuracy on Xs is inaccessible without the Q-Lovc 
relation. 

On a fundamental physics front, independents mea- 
surement of any two members of the I-Love-Q trio would 
allow for model-independent and EoS-independent tests 
of GR. For example, let us assume that one has measured 
the moment of inertia of the primary NS of the double 
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FIG. 3. (Color Online) (Top) I-Love relation in GR (black 
thick solid curve), normalized as in Fig. 1, with A/ and 
^X(tid) giiQ^ji ag (black) dashed lines around a fiducial mea- 
sured value, shown with a cross. AJ is 10% of the fidu- 
cial value, assuming future double binary pulsar observa- 
tions [6, 7]. AA<"'*' is 60% of the fiducial value, assum- 
ing a GW binary NS observation. We also plot the I- 
Love relation in dynamical CS gravity with a CS parameter 
^cons = 1.1 X 10* At* for various EoSs. This test would con- 
strain ^ < ^cons, 6 orders of magnitude more strongly than 
current Solar System bounds [56]. For reference, Af, for the 
Shen EoS is shown on the top axis. (Bottom) Relative frac- 
tional difference of the CS I-Love relations (light solid curves) 
between the realistic EoSs and the n — 1 polytrope. For 
reference, we also show this difference for the GR I-Love rela- 
tions (thick solid curves) . Observe that even in dynamical CS 
gravity, the universality of the I-Love relation seems to hold 
to 0(1) %. 



binary pulsar J0737-3039 to 10% accuracy^ [6, 7]. Let 
us further assume that GW observations have indepen- 
dently measured the NS tidal Love number to roughly 
60% from a detection of an equal-mass NS binary with 
the same NS mass as the primary in J0737-3039. With 
these observations, one can then plot a point in the I-Love 
plane with a measurement error box as shown in the top 
panel of Fig. 3. Such a figure automatically provides a 
consistency (null) test of GR: one can test whether GR 
predicts an I-Love curve that goes through such an error 
box. Moreover, one can also constrain modified gravity 
theories by requiring that the I-Love curves in these the- 
ories pass through this error box. We will show here that 
such a test is even possible when the GW binary system 
has component masses that are up to 10% different from 
the pulsar ones. 



^ Notice that this pulsar has a relatively long spin period for a 
millisecond pulsar, 22.7 ms, and thus, the slow-rotation approx- 
imation is perfectly valid. 



5 



Such a test, of course, is constraining provided modi- 
fied gravity theories predict I-Love-Q relations that are 
not degenerate with the GR ones. Figure 3 shows that 
at least for dynamical Chern-Simons (CS) gravity [57] 
this is not the case. Dynamical CS gravity is a parity- 
violating and quadratic-curvature corrected theory that 
has been weakly constrained in the Solar System through 
Gravity Probe B observations [56] and table-top exper- 
iments [58]. The fiducial I-Love measurements of Fig. 3 
would constrain this theory 6 orders of magnitude more 
strongly than current tests, down to ^ < 1.1 x lO'^Mf, 
where ^^Z'' is the characteristic length scale of the the- 
ory. Observe also that the universality of the I-Love re- 
lation still holds in dynamical CS gravity within several 
% accuracy, although this universality does not hold as 
well as in GR. Perhaps, this is because NSs in dynamical 
CS gravity have scalar hair that depends nontrivially on 
the NS's internal structure [59], and thus, the internal- 
structure effacing in dynamical CS gravity may not be as 
effective. 

The I-Love-Q relations presented here show universal- 
ity within the framework (uniform and slow-rotation and 
small tidal deformations) we work in. Of course, this 
framework is inappropriate to study newly-born NSs, 
which are probably differentially rotating and doing so 
fast. Older NSs that are uniformly rotating usually spin 
slowly, especially those that serve as a source of GWs 
for ground-based detectors, as they will have spun down 
significantly by the time they enter the GW sensitivity 
band [60]. Short-period, millisecond pulsars, on the other 
hand, spin at a non-negligible rate, and thus, a slow- 
rotation expansion may not be suitable. In that case, 
we still expect to find semi-universal I-Love-Q relations 
that although independent of the EoS will depend on the 
spin period. Such an analysis requires a full numerical 
treatment of rapidly rotating NSs [61-64] and is left for 
future work. 

The rest of this paper deals with the details of this 
calculation and it is organized as follows. In Sec. II, 
we explain how the spacetime is decomposed, the ap- 
proximations used and the stress-energy tensor we used 
to describe NSs. In Sec. Ill, we construct non-rotating, 
isolated NS solutions. Taking these solutions as a back- 
ground, in Sec. IV we construct slowly-rotating NS so- 
lutions to linear order in spin and calculate the NS mo- 
ment of inertia. In Sec. V, we construct slowly-rotating 
NS solutions to quadratic order in spin, calculate the NS 
quadrupole moment and define the rotational Love num- 
ber. In Sec. VI, wc define the tidal Love number and 
construct tidally-corrected NS solutions. In Sec. VII, we 
show how the I-Love-Q relations are essentially indepen- 
dent of the NS EoS. We also study explanations for these 
relations, by considering analytic relations when using 
polytropic EoSs in the Newtonian limit. In Sec. VIII, 
we explain possible applications of the I-Love-Q rela- 
tions to observational astrophysics, GWs and fundamen- 
tal physics. In Sec. IX we conclude and point to future 
research. 



All throughout the paper, we follow mostly the con- 
ventions of Misner, Thorne and Wheeler [65]. We use 
the Greek letters (a, /?,••• ) to denote spacetime in- 
dices. The metric is denoted by g^jy and it has signature 
(—,-!-, We use geometric units, with G = 1 = c. 



II. SPACETIME DECOMPOSITION AND 
MATTER REPRESENTATION 

In this paper, we consider uniformly rotating NSs that 
are slightly deformed either due to rotation or tidal fields. 
Such solutions can be numerically constructed pertur- 
batively in a slow-rotation and tidal-deformation expan- 
sion, taking the non-rotating, isolated solution as a back- 
ground. In this section, we explain the metric decompo- 
sition employed here and the stress-energy tensor we will 
use to describe NSs. 



A. Metric Decomposition 

We choose Boyer-Lindquist type coordinates (i, r, 0, ( 
and decompose the metric as 



2e'm2{r)aY2m{e, 



Mr) 



de 

dr^ 



r - 2M{r) 
+7-2 [l + 2^K2{r)aY2^{e,^)\ 

X \d?e + sin^ e [d(j} - €[VL^ - uii{r)P[{cosd)]dtf'^ 

+0{e^), (1) 
where M{r) is defined by 



M{r) 



1 - e-^W 



(2) 



Pt{cos9) is the £-th order Legendre polynomial, and 
^£m(^, 0) is the spherical harmonic function. The quan- 
tity e here is a book-keeping parameter that we will later 
set to unity and we only introduce to remind ourselves 
of the order of the approximation. Terms linear in e 
are induced only by linear-order in rotation effects, while 
tidal-deformation effects enter at O(e^). We will work 
here to quadratic order in e. 

A slow-rotation expansion is quite appropriate to 
model old NSs. Recycled millisecond pulsars, the fastest 
NSs observed to date, have angular velocities in the 
kHz, but this number is small relative to the NS mass, 
i.e. M,il* < 0.01, where f2* is the NS angular velocity. 
For the fastest millisecond pulsar J19394-2134 [66], with 
period 1.5 ms, the dimensionless spin parameter, defined 
via X = S/Ml = m^/M^, is still small x ^ 0-3, using a 
Newtonian expression for the moment of inertia. Thus, a 
slow-rotation expansion is well-justified, especially when 
carried out to second order. This approximation, how- 
ever, would break down if considering newly-born NSs, 
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which are likely to be differentially rotating, much hotter 
and with much larger magnetic fields. Notice also that 
the NSs that will source GWs in the band of ground- 
based detectors are expected to have significantly smaller 
spins than that. This is because NSs spin-down [60] as 
they inspiral and ground-based detectors will only be sen- 
sitive to the last 17 minutes of the orbit before coales- 
cence. 

The free functions in our metric decomposition are v 
and A at 0(e°), Wi at 0{e) and h2, K2 and TO2 at 0{e'^). 
The leading-order correction due to slow rotation enters 
at 0{e), while that due to tidal deformations enters at 
O(e^). For the former, we restrict ourselves to axisym- 
metric perturbations; at 0{e) only the {£,m) = (1,0) 
mode survives, while at O(e^) only the {i,m) = (0,0) 
and {£,m) = (2,0) modes survive. For the latter, we are 
only interested in the spin and tidal, quadrupolar defor- 
mations, and thus we only keep (. = 2 modes in Eq. (1), 
but allow for all m modes. Henceforth, we set the con- 
stant a = 2-\/7r/5 so that aY2o{9, 4>) = P2{cos9). 

As pointed out by Hartle [37], one needs to be care- 
ful about choosing coordinates when deriving and solv- 
ing perturbed equations. A perturbative analysis is valid 
only if perturbed quantities are much smaller than the 
unperturbed one. If one were to carry out calculations 
in {t, r, 9, </>) coordinates, such conditions would be vio- 
lated in certain situations. For example, in the region 
of spacetime outside the unperturbed star but inside the 
perturbed star, the ratio of the perturbed pressure (or 
density) to that of the unperturbed pressure (or density) 
diverges, which violates our perturbative treatment. 

In order to overcome this problem, we transform the 
radial coordinate via [37] 

r{R, 9)=R + e^URhY2ra{9, </)) + ©(e^) , (3) 

where ^2 (R) is such that 

p[r{R,e,c^)]=p{R) = p^°HR). (4) 

In other word, the new radial coordinate R is chosen 
such that p[r{R^ 6, (j))] is identical to the unperturbed den- 
sity p*-"-* (r) . By construction, the density and pressure in 
these new coordinates contain only the unperturbed con- 
tributions. Notice that ^2i^2?n is well-defined only inside 
the star and we take it to be constant outside. This 
means that the exterior metric in (<, r, 0, (p) coordinates 
can be obtained simply by replacing i? ^ r in the exte- 
rior metric in {t, R, 9, 4>) coordinates. 

The transformed metric in (t, R, 9, (j)) coordinates can 
be found in [67] for the axisymmetric case. Henceforth, 
we will relabel the metric coefficients via 

iy{R) = D{r) = v{R + e^S,2aY2ra) , 
X{R) = A(r) = \{R + e^C2aY2^m) , 
wi {R) = wi (r) =u}i{R + €^£,2aY2„,) , 
h2{R) = hir) = h2{R + <?iiaY2^) , 
7712 (-R) = 7712 (r) = m2(-R + e^^2aY2m) , 



K2{R) = K2{r) = K2{R + e'6«>2m) , 

M{R) = M{r) = M{R + e^£,2aY2^) . (5) 

B. Matter Representation 

We here consider NSs that are uniformly rotating, and 
thus, we model them as a perfect fiuid. Uniform rotation 
should be a reasonable approximation unless one consid- 
ers newly-born NSs. The stress energy-momentum tensor 
of the matter field T^'"' is then given by 

r/r' = {P + P)uf,u„ + p g^^ , (6) 

where the four-velocity is given by 

li^ (w°,0,0,el7,u"), (7) 

and fl* is the constant angular velocity of the NS. By 
using the normalization condition w^w^ = —1, we obtain 
the time component of the four- velocity vP as 

w° = e-''/'^ + e''—^[ujlP[^R^ Bin^ 

-e''{2h2 + y'i2)aY2^]+0{e'). (8) 

We here consider 4 realistic EoSs: APR [38], SLy [39, 
68], Lattimer-Swesty with nuclear incompressibility of 
220MeV (LS220) [40, 69] and Shen [41, 42, 69], the lat- 
ter two with temperature of O.OlMeV and an electron 
fraction of 30%. All of the EoS described above are "re- 
alistic" in that they allow NSs with masses larger than 
1.93Mo, the lower bound of the recently found massive 
pulsar J1614-2230 [70]. For comparison purposes, we also 
consider polytropic EoSs, i.e. EoSs of the form 

p = A>i+i/" , (9) 

where K is an amplitude constant and n is the constant 
polytropic index. One can approximate the NS EoS with 
polytropes in the range n « 0.5 — 1 [17, 71]. No single 
polytrope, however, is believed to be an accurate repre- 
sentation of a realistic EoS. 

The APR EoS is constructed by using the variational 
chain summation methods, which is expected to include 
all leading many-body correlation effects. The APR EoS 
uses Hamiltonians that include a three-nucleon interac- 
tion, which predicts that a transition exists from NS mat- 
ter to a phase with neutral pion condensation at a baryon 
number density of ~ 0.2 fm~^. The SLy EoS is calculated 
from a non-rclativistic mean field theory approach, with 
a new set of Skyrme-type effective nucleon-nucleon inter- 
actions, suitable for describing very neutron rich matter. 
Unlike the APR EoS that describes only the NS's liquid 
core, the SLy EoS is a "unified EoS" in the sense that 
it is supposed to describe also the NS crust. The LS220 
EoS is constructed from a finite-temperature compress- 
ible liquid-droplet model with a Skyrme nuclear force. 
Such an EoS is derived within the single heavy nucleus 
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approximation and the assumption of nuclear statisti- 
cal equilibrium. The Shen EoS uses a relativistic mean- 
field theory model and assumes nuclear statistical equilib- 
rium. Nuclear incompressibility of the Shen EoS occurs 
at 281MeV. 



III. SLOWLY ROTATING, ISOLATED 
NEUTRON STARS: 0(x°) 

In this section, we construct non-rotating, isolated NS 
solutions, which will later be used as background solu- 
tions to construct slowly-rotating, tidally-deformed NS 
solutions in Sees. IV, V and VI. 



A. Einstein Equations and Exterior Solutions 

The {t, t) and {R, R) components of the Einstein Equa- 
tions yield 



df _ AttR^p + M 
dR " R{R - 2M) ' 



(10) 
(11) 



respectively. Combining the i?-componcnt of the equa- 
tion of motion VT^"^' = and Eq. (11), one obtains the 
Tolman-Oppenheimer-VolkofF (TOV) equation: 



dp _ {4ttR^p + M){p + p) 
dR ~ R{R - 2M) 



(12) 



Equations (10), (11) and (12) together with the equation 
of state p = p{p) close the system of differential equa- 
tions. 

The exterior solutions to the above equations can be 
obtained by setting p = Q = p. One finds [37] 



^ext(^^^) = -A<=''*(i?) = In 1 



R 



(13) 



We use the superscripts "ext" to refer to exterior quan- 
tities. 



B. Interior Solutions 

First, we solve Eqs. (10) and (12) together with the 
equation of state with initial conditions 



p(r,) = + 0{rl) , 
p{r,:)^p, + 0{rl), 

4-71 

M{r,) = -py^+0{r!), 



(14) 
(15) 

(16) 



where pc and p^ are the central density and pressure re- 
spectively, and corresponds to the core radius which 
we take to be = 100 cm <c 7?.*. We have checked 




FIG. 4. (Color Online) The NS mass-radius relation for a few 
realistic EoSs. The black, horizontal, dashed line at 1.93M0 
corresponds to the lower bound on the mass of the recently 
found massive pulsar J1614-2230 [70]. Observe that all real- 
istic EoSs lead to mass-radius curves that exceed this lower 
bound. 



that all of our results are independent of the choice of 
provided this is a very small number relative to the NS 
radius. We solve Eqs. (10) and (12) outwards from r = 
until p vanishes. The NS radius 7?.* and the NS mass M^, 
are then defined by p{TZ^) = and M* = MiJZ^) re- 
spectively. For later convenience, we introduce the NS 
compactness 



C 



7^* 



(17) 



Notice that The central pressure Pc = p{Pc) is determined 
from the EoS, once pc is chosen. The central density pc 
is then a free parameter that effectively determines the 
mass and radius of the NS. 

With this solutions, we can then solve Eq. (11). One 
approach is to use the boundary condition [see Eq. (13)] 



1 



2Ah 

lz7 



(18) 



at the NS surface as an initial condition and then in- 
tegrate inwards toward the core. Another approach is 
to use the fact that Eq. (11) is shift invariant, as done 
e.g. in [67]. All throughout this paper, numerical solu- 
tions to the initial value problem are obtained with an 
adaptive 4th-order Runge-Kutta method [72]. 

Figure 4 shows the mass-radius relation for various 
EoSs. We have checked that the mass-radius relation for 
the SLy EoS agrees with that shown in [68]. As antici- 
pated, all EoSs lead to NSs with maximum mass larger 
than 1.93 Mq (the black dashed horizontal line), which is 
the lower bound for the mass of J1614-2230 [70]. We do 
not show the mass-radius relation for the polytropic EoSs 
because the I-Love-Q relations that we present in Sec. VII 
only depend on the NS compactness and do not depend 
on the mass-radius relation. Figure 5 shows the interior 
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FIG. 5. (Color Online) The NS density (top) and pressure 
(bottom) profile for various EoSs as functions of the radial 
coordinate. We fix C = 0.17 which corresponds to M, = 
1.4M0 for the APR EoS. We do not show profiles of n = 2.5 
and n = 3 polytropic EoSs since the maximum compactness 
with such EoSs is smaller than 0.17. 



profile of the NS density (top) and pressure (bottom) as 
functions of the radial coordinate for a compactness of 
C = 0.17, which corresponds to a NS with M, = IAMq 
and « 12.1 km for the APR EoS. 



IV. SLOWLY ROTATING, ISOLATED 
NEUTRON STARS: O(x^) 

Let us now focus on constructing slowly-rotating, iso- 
lated NS solutions. In this section, we only consider ax- 
isymmetric perturbations and construct NSs to linear or- 
der in spin. We will first discuss the differential equation 
that needs to be solved, and then we will solve them in 
the exterior region modulo an integration constant. After 
this, we discuss the asymptotic behavior of the solution at 
the NS center, which can then be used as an initial condi- 
tion to solve the equations in the interior region. Finally, 
we determine the integration constant by matching the 
interior and exterior solutions at the NS surface. 



A. Einstein Equations and Exterior Solutions 

At linear order in e, the only non- vanishing component 
of the Einstein Equations is the {t, (f) one: 



dR? ^ R IlR 



- 167r(p + p)e^wi = 0. 



(19) 



Solving this equation in the exterior (i.e. setting p 
p), one finds [37] 



cxt 



2S_ 

R3 



21 
R^ 



where we have defined the moment of inertia by 



(20) 



(21) 



This quantity characterizes how fast a body can spin 
given a fixed spin angular momentum S. Notice that 
the exterior solution depends on two constants D,^ and 
S. The former must be specified a priori, just like pc, 
and it describes how fast the NS is rotating. The latter 
is determined by matching this exterior solution to an 
interior solution at the NS surface. 



B. Interior Solutions 

Before we can solve for the interior solution, we 
first need initial conditions at the NS center. Taylor- 
expanding Eq. (19) about the NS center, we find that 
the interior solution must asymptotically behave as 



t^l(i?) =Wc+— (Pc+PcRi?'+0(i?') (R 



0+). 
(22) 

This solution contains a single constant, lUc, because we 
have eliminated another constant by requiring regularity 
of the solution at the NS center. The constant ojc deter- 
mines the NS spin angular momentum S, or equivalently, 
the NS moment of inertia /; in particular, / increases as 
LOc increases. 

We numerically solve Eq. (19) with the initial condi- 
tion in Eq. (22) via an adaptive 4th-order Runge-Kutta 
method [72]. In solving this equation, one can take ad- 
vantage of its homogeneity, its scale- invariance, as done 
e.g. in [67, 73]. Once the interior solution has been found, 
we match it to the exterior one in Eq. (20) at the NS sur- 
face R = TZt, . The matching ensures that the solution is 
continuous and differentiable at the NS surface: 



int 



(7^.) 



/ int 



(7^.) = ^^°-'(7^*), (23) 



where we use the superscript "int" to refer to interior 
quantities. Through these conditions, we determine S 
(or equivalently /) and Wc as a function of 17*. In prac- 
tice, due to the scale invariance of Eq. (19), the exterior 
solution can be divided by fi* and thus it only depends 
on the single constant /. Similarly, the interior solution 
can be obtained for w™*/Sl* as a function of a single con- 
stant uic = Wc/O*. Therefore, the conditions in Eq. (23) 
uniquely determine / and ujc- This then determines the 
full solution, and thus also S, up to the overall constant 
of proportionality Q^. 
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FIG. 6. (Color Online) Dimensionless moment of inertia J, defined in Eq. (26), as functions of M, (left) and C (right) for 
various EoSs. The horizontal dashed lines at / = 21.1 correspond to M« = IMq for the APR EoS; NSs below this line have 
higher M* and C. The solid cross indicates the value of I for a BH. Observe that the / curves for realistic EoSs approach each 
other as C increases, and moreover they approach the BH limit 7bh = 4 as C — >• 0.5. 



The moment of inertia can be expressed entirely as 
a function of the interior solution. From Eqs. (10)- 
(19), (20) and (23), / takes the form [37, 74] 



cross in Fig. 6. Of course, none of the NS sequences con- 
sidered here will ever lead to a BH solution for any finite 
choice of central density. 



3 Jo R- 2M{R) 



(24) 



In the Newtonian limit (superscript "N"), Eq. (24) re- 
duces to [37] 



Stt 



r = — / R^p{R)dR, 



(25) 



For later convenience, we define the dimensionless mo- 
ment of inertia / 



V. SLOWLY ROTATING, ISOLATED NEUTRON 
STARS: O(x^) 

Let us now look at slowly-rotating NS solutions at 
quadratic order in spin. Following Sec. IV, we first dis- 
cuss the differential equations that describe the solution 
and then we solve them in the exterior region. We then 
discuss the asymptotic behaviors of the solutions at the 
NS center, obtain the interior solutions numerically, and 
match it to the exterior solution at the NS surface. 



MS 



(26) 



A. Einstein Equations and Exterior Solutions 



Figure 6 shows / as a function of the NS mass 
and compactness C. We have verified that the moment 
of inertia / obtained here agrees exactly with previous 
results in the literature [75]. Observe that the different 
/ curves for realistic EoSs approach each other as C in- 
creases. Moreover, observe that all these curves approach 
the value of / for a BH as C — > 0.5, shown with a solid 

I 



At quadratic order in spin, the ^-component of the 
equation of motion V^T^'g' = 0, valid only inside the 
star, yields 



R:^e~^{3h2 + e-"'R'^ujl) 
~ 3(A/ + AnpR^) 



The 



and (R, R) components of the Einstein Equations give respectively. 



7712 

dK2 
dR 
dh2 

Ir 



-Re-^h2 + -R'^e-^''+^^ 
6 



Re 



dh2 R-iM - 

Hr ^ W 

R- M + AirpR 



^T^PR" A, 

e 



dR 

R-M 



AnpR' 2X 



i?3 



3 ^dK2 



R 



dR 



3-47r(p-Hp)i?2 
R 



m2 ; 

R 



C2 + 



SirpR' 



i?2 



12' 



(27) 



(28) 
(29) 



duji 

Ir 



10 



47r(p + p)i ?V 
3i? 



By imposing asymptotic flatness at spatial infinity, one finds the exterior solutions [37] 



1 



1 

1 



1 
1 



2ilL 



.2 , 



- 2M, 



1 - 7- 



R 



3AR 



10- 



5^ 



A/* 2 Ml 
^ ~r" ^ 3"r2" 



R 
2MI 



1 



A/, 



A/, 4 A/2 
^-^ir + 3# 



M 
2A/2 

3 i? 3 



ln/(i?) 



2A/* 



/(i?)^ln/(i?) 



(30) 



(31) 

(32) 
(33) 



with f{R) = 1 - 2Af*/i?, and the associated Leg- 
endre functions of the second kind and A an integration 
constant that is to be determined by matching with the 
interior solution at the NS surface. 

The spin- induced quadrupole moment (5^''°*) can be 
read off from the coefficient of the P2{cos9)/ R^ term in 
the Newtonian potential [76]: 



^(rot) 



-AMt 



(34) 



Notice that the quadrupole moment depends both on 
the integration constant S and the integration constant 
A, determined after matching the interior and exte- 
rior linear- and quadratic-order in spin solutions at the 
NS surface. The quadrupole moment represents the 
quadrupolar deformation of a body away from sphericity, 
with (3^'°*-* < corresponding to an oblate deformation. 
Notice also that the first term of Eq. (34) is identical to 
the relation one obtains for BH, which means that ^ — > 
in the GR test-particle limit. 



B. Interior Solutions 

Let us begin by Taylor-expanding Eqs. (27)-(30) about 
the NS center and solving the expanded equations to ob- 
tain 



h2{R) = BR^+OiR^), 
K2{R) = -BR^ + 0{R'^), 
m2{R) = --BR^ + 0{R^) , 
35 + e-^-uj^ 



47r(pc + 3pc) 



R + 0{R'^), {R^Q-^ 



(35) 
(36) 
(37) 



(38) 



where i? is a constant that determines the NS quadrupole 
moment. As before, the constant is defined as = 
i/(re). 



We numerically solve the evolution Eqs. (29) and (30) 
with the initial conditions of Eqs. (35) and (36), using 
an adaptive 4th-order Runge-Kutta algorithm [72]. As 
before, when solving these equations we must impose the 
following boundary conditions, such that /12 and K2 are 
continuous at the NS surface: 



;^■"'(7^.) = /^r(7^*), 



(7^. 



= if2"^*(7^*). (39) 



These matching conditions determine the constants A in 
Eqs. (31) and (32) and B in Eqs. (35) and (36). 

In practice, we follow [37, 67] and first solve the interior 
solution as a sum of a particular solution, with some 
test-value for B, and the product of an undetermined 
constant and the homogeneous solution. We then fix this 
undetermined constant, together with A, by requiring 
that the interior and exterior solutions match at the NS 
surface. We have checked the results obtained through 
this method by solving the equations using the Riccati 
method [77-79]. 

Figure 7 shows the dimensionless rotationally-induced 
quadrupole moment Q as functions of A/* and C, where 
Q is defined by 



Q 



(40) 



where we recall that the dimensionless spin parameter 
X is defined by x = S/M^. This Q is the same as the 
dimensionless quadrupole moment a in [80]. As in the / 
case, the Q curves for realistic EoSs approach each other 
as C increases. Moreover, these curves also approach the 
Q value for a BH as the compactness approaches 0.5. As 
before, however, the NS sequence does not go to a BH 
solution for any finite choice of central density. 



C. Rotational Love Number 

With the quadratic isolated NS solutions at hand, we 
can now introduce the rotational Love number [81]. In 
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FIG. 7. (Color Online) Dimensionless quadrupole moment Q, defined in Eq. (40), as functions of M, (left) and C (right) for 
various EoSs. The horizontal dashed lines at Q = 9.89 correspond to a NS with M* = IM©; curves below this line have higher 
M« and C. Observe that the Q curves for realistic EoSs approach each other as C increases, and moreover, they approach the 
BH limit Qbh = 1 as C -S- 0.5. 



general, Love numbers represent the deformability of a 
NS away from sphericity. The rotational Love number, 
in particular, refers to the deformability of a NS due to 
its spin. 

Love numbers are defined in a buffer zone, the region 
TZ ^ R ^ ^Z^,, where TZ is the radius of curvature of 
the source of the perturbation. For example, the {t, t) 
component of the metric can be expanded in the buffer 
zone as [8, 81-85] 



gtt 



47r 
'l5 



5 i?3 



i?3 

7^3 



(41) 



The quantity is related to the trace of the 

rotationally-induced, electric, quadrupole tidal tensor, 
i.e. the quadrupolar contribution of the centrifugal po- 
tential. In the Newtonian limit, this quantity reduces to 
£.(rot) ^ ^2 yg^j^i^ Y2m{^) are the € = 2 spheri- 

cal harmonics in the Cl direction, where h is the principal 
axis of the perturbation, which in this case corresponds 
to the unit vector of the spin angular momentum S. 

The (. = 2 rotational Love number A(''°') is defined 
by [81, 86] 



Q(rot) 
" £:(rot) 



2 (rot) 



(42) 



where the second equality uses the Newtonian expression 
for S^™*). As defined here. A'™'' has unit of (mass)^ or 



(length)^ (recall that we use geometric units throughout 
this paper, where c ^ 1 = G), and thus, there are 2 
natural ways of normalizing it [81, 86]; 



^(rot) _ 



3 A(™*) 

2 ni 



\(rot) 9 , 

Af 5 3 2 



(43) 
(44) 



By using Eqs. (24), (26) and (40), one can rewrite A^™*) 
as 



72, 



In this paper, we refer to ^2'^°''' as the (. = 2 rotational 
apsidal constant, while we refer to A^™') as the ^ = 2 
dimensionless rotational Love number. 



(45) 



VI. TIDALLY-DEFORMED NS SOLUTIONS 

Up until now we have concentrated on isolated NSs 
in the slow-rotation approximation. We will now switch 
gears and consider NSs in a binary system. We focus on 
one of the binary components, the primary, and study 
how it is tidally deformed by its companion, assuming 
the primary is not spinning. One can construct tidally- 
deformed NS solutions in a manner similar to the con- 
struction of slowly-rotating solutions. In both cases, the 
deformation (either due to rotation or tidal effects) is 
treated as a small deformation away from sphericity. 



A. Einstein Equations and Exterior Solutions 

The leading-order effect of tidal perturbations enters 
at 0{e'^). This is because this effect is generated by an 
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electric tidal perturbation, which must be of even par- 
ity. Moreover, in this section we are interested in non- 
rotating tidally deformed NSs, so we can set wi = in 
Eqs. (28)-(30). By eliminating m2 and K2 from these 3 
equations, one obtains a master equation for /12 [8]: 



2M 

— + 47ri?(p - p) 




9p- 




The observable related to tidally-deformed NS will 
eventually be a tidal Love number, and thus, we will 
need to asymptotically expand the exterior solution in 
the buffer zone. This time, however, the radius of cur- 
vature that defines the buffer zone is related to the tidal 
field generated by the companion. This radius is ap- 
proximately equal to the orbital separation of the binary. 
Therefore, when solving Eq. (46) in the exterior region, 
one cannot impose asymptotic flatness to eliminate one 
of the constants of integration. Keeping this in mind, the 
solution to the above equation is [8] 



hT' = ci 



1 



R 



2Ah 
R 



R^{R - 2Ah 



-3 In 



R 



R - 2A'L 



C2 



R 



1 - 



2AL 



R 



(47) 



where ci and C2 are integration constants. 



B. Interior Solutions and the Tidal Love Number 

The interior solution to Eq. (46) can be obtained by 
solving this equation numerically with the initial condi- 
tion in Eq. (35) and its derivative. We obtain this nu- 
merical solution in the same way as we obtained /i2 for 
slowly-rotating NSs. As before, the interior solution will 
depend on the integration constant B, which, in princi- 
ple, is determined by matching this solution to the exte- 
rior solution in Eq. (47) at the NS surface: 

h^'HT^*) = /^r(7^*), h'2"''in,) = /^'2'="*(7^,) . (48) 

Notice that by using Eq. (48), wc can re-express ci and 
C2 in terms of h2(TZ*), /i2(^*) ^^'^ the NS compactness 
C. 

With the interior solution in hand, let us now define 
the tidal Love number. As in the case of the rotational 
Love number, the tidal one characterizes the deforma- 
bility of a NS away from sphericity, but this time due 
to the presence of a tidal field induced by a companion. 
In the buffer zone, the (t, t) component of the metric 



takes the form of Eq. (41), but with Q^™*) ^ Q'^^''^^ and 
£(rot) ^ £.(tid)^ ^j^g^g g(tid) £.(tid) correspond to the 

tidally-induced quadrupole moment and the tidal poten- 
tial, induced by the companion. We then define the tidal 
Love number A^'"^^ by 



^(tid) 



g(tid) 



and its dimensionless versions ^ ' and A("^) by 
zJtid) _ 3 A^"'^) 



2 n 

A(«d) 

M5 ■ 



(tid)^-5 



(49) 

(50) 
(51) 



Following [17], we here refer to A^*"^' as the tidal Love 
number''* and wc refer to fcj'"^'' and A^*"^^ as the tidal 
apsidal constant and dimensionless tidal Love number, 
respectively. 

The prescription of the tidal Love number is completed 
by finding the value of f which is determined by the 
asymptotic behavior of /12 in the buffer zone. Taylor- 
expanding this quantity in the buffer zone, one finds [8] 



16 



C2 



O 



Mf R 



(52) 



As shown in Eq. (41), the term in the asymptotic expan- 
sion of gtt (or /12) in the buffer zone that is proportional 
to R~^ gives us the tidal quadrupole moment, while the 
term proportional to R^ gives us the tidally-induced elec- 
tric quadrupole tidal tensor. Thus, we find that ci is 
related to Q^*'"^), while C2 is related to £^'^''^\ 

The tidal apsidal constant can then be found by taking 
the ratio of ci and C2 [8] : 



Z,(tid) 



^C^El = |c'5(l _ 2Cf[2 + 2C{v - 1) 



V] 



5 C2 5 
X {2C[6 - 3y + 3C(5j/ - 8)] 

-|-4C^[13 - lly + C{?,y - 2) + 2C^(1 + y)] 

-f3(l - 2Cf[2 -y + 2Ciy - 1)] ln(l - 2C)}"' , 

(53) 

with y = Ti^:h'2{TZ*) / h2(Ti-*) ■ In the second equality, we 
have rewritten ci,2 in terms of /12, its derivative and the 
NS compactness. 

We see then that the tidal apsidal constant only de- 
pends on y, which simplifies the way one must solve 
Eq. (46). First, wc notice that Eq. (46) is a homoge- 
neous equation for h2, and thus, the integration constant 
B in Eq. (35) only changes the solution h2 by a constant 



^ In some references, A'*"^) is called the tidal deformability and 



the word "tidal Love number" is reserved for k. 
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FIG. 8. (Color Online) Dimensionless tidal Love number A*-*"^^ as functions of M* (left) and C (right) for various EoSs. The 
horizontal dashed lines at A''"*' = 2.66 x 10^ corresponds to a star with M* = 1 Mq; the region below this line corresponds to 
stars with larger mass and compactness. Observe that the A^*"^^ curves for realistic EoSs approach each other as C increases, 



and moreover, approaches the BH limit A 



(tid) 



as C ^ 0.5. 



factor. Since y cx /i2/ft.2 does not depend on this overall 
factor, it suffices to solve Eq. (46) with an arbitrary test 
value for B, if one is only interested in the tidal apsidal 
constant. We have calculated the tidal apsidal constant, 
as well as the tidal Love number for a sequence of stars 
with varying M* and C . We have found that our results 
agree exactly with Figs. 1 and 2 of [19]. 

Figure 8 shows the dimensionless Love number A'-'"^-' 
as functions of Af, and C. Observe that the A^'"^^ curves 
for realistic EoSs approach each other as C increases, and 
moreover, they approach the BH limit as C — J- 0.5. Once 
more, as before, the BH limit cannot be taken from the 
sequence of NS considered, as there is no finite central 
density that would lead to BH formation. 

VII. I-LOVE-Q RELATIONS 

Now that the moment of inertia, quadrupole moment 
and Love numbers have been calculated, let us present 
the universal I-Love-Q relations. We first show numerical 
results and a fitting curve through these. Then, we obtain 
analytic I-Love-Q relations for the n — Q and 1 polytropic 
EoSs in the Newtonian limit. 



A. Numerical Results 

Figure 9 shows universal relations between dimension- 
less quantities, /, Q, A'*"^^ and A'™*' for various EoSs. 
Notice that these (barred) 4 quantities are independent of 
the NS spin, to second order in the slow-rotation approx- 
imation. The parameter varied along each curve is the 
NS central density, or equivalently the NS compactness. 
Therefore, for the polytropic EoS the I-Love-Q relations 
are independent of the polytropic amplitude coefficient 



K in Eq. (9). The bottom part of each panel shows the 
relative fractional difference between each of the curves 
and curve corresponding to the n — 1 polytropic EoS. For 
reference, the top axes show the NS mass with the APR 
EoS. The vertical dashed lines correspond to A/* — IMq 
for the APR EoS and points to the left of these lines 
correspond to more massive NSs with higher compact- 
ness. Observe that, for realistic EoSs with Af* > IMq, 
the fractional relative differences are 0(1)%. Observe 
also that the polytropic I-Love-Q relations deviate from 
those with realistic EoSs as one increases n, ie. as the 
NS becomes more centrally-concentrated. We see this 
as evidence that the I-Love-Q trio is most sensitive to 
the NS outer layers, where realistic EoSs mostly agree 
with each other. Curiously, the fractional relative differ- 
ence between the n = 1 and n = polytrope (constant 
density NS star) is also of 0(1)%, but this case will be 
studied analytically in Sec. VII B. 

Observe that, in general, the dependence of the I-Love- 
Q relations on any EoS becomes weaker as the NS com- 
pactness C increases, ie. from right to left in any of the 
panels of Fig. 9. This may be evidence that at least part 
of the universality observed is due to the NS sequence 
approaching a BH as C — >■ 0.5, where the latter has no 
internal-structure dependence by the no-hair theorems. 
As C — )- 0.5, the asymptotic values of /, A*-'"^-* and Q for 
a BH_are 7^4 [87], A^"'^) ^ [10] [see also Eq. (53)] 
and Q — )• 1 [49, 50, 80] respectively. In each panel of 
Fig. 9, wc show the BH limit of either / = 4, Q = 1 
or A(™') ^ 16 as a horizontal dashed line. Observe that 
the I-Love-Q relations asymptote to such BH values as C 
increases. However, one cannot quite reach this limit, as 
one can never construct a BH solution by increasing the 
central density of a NS solution by a finite amount. We 
think that this is why the relative fractional differences 
shown in the bottom panels of Fig. 9 do not decrease to 
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FIG. 9. (Color Online) I-Love (top left), Q-Love (top right), I-Q (bottom left) and Love-Love (bottom right) relations for 
various EoSs. Black thin solid lines represent the Newtonian limit for the polytropic EoSs with n = and n = 1. The 
horizontal lines at / = 4, Q = 1, A*-*"*' = and A*-"^"*' = 16 correspond to the (non-rotating) BH limiting values. Observe 
that as one increases the NS compactness (toward the left of each panel), the I-Love-Q relations approach the BH limit. The 
(barred) quantities I, Q, A*-*"*' and A'"^"'-* do not depend on the NS spin to second-order in the slow-rotation approximation. 
The parameter varied along each curve is the NS central density, or equivalently the NS compactness, both increasing to the 
left of the plots. For reference, the vertical dashed lines correspond to M = IMq for the APR EoS. The top axis of each panel 
shows the corresponding NS mass for the APR EoS. The bottom part of each panel shows the relative fractional differences 
between the relations, using the n = 1 polytropic curve as a reference. 



zero as one increases the central density. 

Unlike the other relations, the A'™*^-A(*"^^ relation 
(bottom, right panel of Fig. 9) depends very weakly on 
the EoS even when the NS compactness is relatively large 
(as one approaches the Newtonian limit). In fact, one 
can exactly show that the relation A*-'"*^ = A*-*"^^ holds 
exactly in the Newtonian limit for any EoS, as we will 
discuss in Sec. VII B 1. 

Given the universality of the I-Love-Q relations, one 



can fit them all with a single curve: 

Inj/i ^ a, + hi Inxi + c,{lnxif + di(lnxi)^ -|- et{\nxi)'^ , 

(54) 

where the coefficients are summarized in Table I. Fig- 
ures 1 and 10 show the fitting curves for the I-Love, Q- 
Love and I-Q relations, together with the relative frac- 
tional difference between the fitting curves and all other 
EoS curves. For the polytropic EoSs, we do not show the 
results when n — 0, 2, 2.5 and 3 since such EoSs do not 
model NSs well; instead, we add results when n = 0.6 
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Xi 


bi Ci 

\J/1 


di 




I 




1.47 0.0817 0.0149 


2.87 X 10"* 


-3.64 X lO"-"^ 


I 


Q 


1.35 0.697 -0.143 


9.94 X 10"^ 


-1.24 X 10"^ 


Q 


^(tid) 


0.194 0.0936 0.0474 


-4.21 X 10"^ 


1.23 X 10"* 



TABLE I. Estimated numerical coefficients for the fitting for- 
mula of the I- Love, I-Q and Q-Love relations given in Eq. (54). 

M (APR) [MJ 

2.18 1.87 1.47 1.23 1.0 0.67 



X APR 
SLy 
O LS220 
□ Shen 

Poly trope (n=l) 
Polytrope (n=0, 
• Polytrope (n=0.6) 
— fit 




FIG. 10. (Color Online) (Top) Fitting curve (solid curve) 
and numerical results (points) of the I-Q relation with various 
EoSs. (Bottom) Fractional errors between the fitting curve 
and numerical results. 



and 0.8. As one can see, the fitting curves are accurate 
to within 0(1)% accuracy. 



B. Analytical Explanations 

The universal I-Love-Q relations presented in the pre- 
vious subsection are quite intriguing, and thus, they beg 
for an analytic explanation. We will attempt one here, by 
investigating these relations for certain EoSs that allow 
for an analytical treatment. In particular, we will study 
the n = and n = 1 polytropic EoSs in the Newtonian 
limit, for which the moment of inertia, the quadrupole 
moment and the Love numbers can be computed fully 
analytically. We can then derive the I-Love-Q relations 
analytically as well to try to obtain an analytical expla- 
nation for these relations. 



ishes, ie. Vp x Vvf = 0. Thus, surfaces of constant p and 
^' coincide. One can express such surfaces in terms of a 
radial parameter a as [81] 



r(a, 6,(1)) = a 



(55) 



where fe is the dimensionless distortion function of the 
constant a. This function is related to the £ = 2 tidal 
apsidal constant ^2*"^'''" by 



^(tid),N ^ 3 - ri2{a*) 
2 2[2-|- 772(a*)] ' 

where a = a, denotes the surface of the star and 



'72(a) 



din a 



(56) 



(57) 



This function can be obtained by solving the Clairaut- 
Radau equation [81, 88] 



l)-6 = 0. 



(58) 



with the boundary condition 772(0) = 0, where I?(a) = 
p{a)/p{a) with p representing the mean density of the 
star. 

In the Newtonian limit, both rotational and tidal ap- 
sidal constants can be calculated from Eq. (58) with the 
same boundary condition [81], and hence. 



^(tid),N ^(rot),N 



(59) 



This shows that the rotationally-induced and tidally- 
induced NS deformabilities are exactly the same in the 
Newtonian limit. In GR, non- linear effects modify this 
relation and break the equality. Equation (59) is shown 
as a black, thin, solid line in the bottom, right panel of 
Fig. 9. Notice that all the curves approach this New- 
tonian result as one decreases the compactness, as ex- 
pected. 

One can also calculate the i = 2 tidal apsidal constant 
fc^*"^^''' by taking the Newtonian hmit {R > M{R), p > 
p) of Eq. (53) [8], where Eq. (46) in the Newtonian limit 
is given by 



2 dh2 
R'dR 



A'np 



dp 
dp 



6 

i?2 



/l2 = . 



(60) 



1. Newtonian A for Generic EoSs 



2. Polytropic I-Love-Q relations: n = 



In Newtonian theory, the curl of the equation of hy- 
drostatic equilibrium, Vp = pV^ where 4" is the total 
gravitational potential, with the Newtonian force van- 



Now, let us investigate the I-Love-Q relations in the 
Newtonian limit for specific EoSs. First, we focus on the 
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polytropic EoS with n = 0, which corresponds to the 
incompressible EoS witli 



Pc 



3 A'h 
47r m 



(61) 



By substituting Eq. (61) into Eq. (25), one obtains 
2 o 2 1 



P = -M 7? 

5 * *' 5C2 



(62) 



Not surprisingly, this is the Newtonian moment of inertia 
for a sphere of constant density. 

Next, we solve Eq. (60) to obtain A^"'^)'". As pointed 
out in [9], one must be careful when solving Eq. (60) 
for the incompressible EoS. This is because p can be ex- 
pressed as a step-function with a discontinuity at the 
NS surface. Thus, dp/ dp in Eq. (60) gives a delta- 
fimction centered at the NS surface. To be more spe- 
cific, by using p ~ —pcS{R — i?*), p(7?.*) = 0, = 
{4:Tt/3)TZIpc and the hydrostatic equilibrium equation 
[or the Newtonian limit of Eq. (12)] at the NS sur- 
face, dp{TZ^,)/dR — — A/*p(72.H.)/7?.J, the coefficient of 
/i2 in Eq. (60) that is proportional to dp/dR becomes 
A-Kpdp/dR = {3/TZ^.)d{R — 7?.*) near the surface. By tak- 
ing such term into account, one can obtain the correct 
^(tici)_N gj,g^ solving Eq. (60) with dp/ dp ~ and then 
shifting by -3 [9]. By taking all of this into account, 

one obtains = -1, and thus, 4"'*^'" =3/4 [9]. This 
agrees with the classic result in [88]. By using Eq. (51), 
A(tid),N becomes 



^(tid), 



1 1 

2C^ 



(63) 



Let us now move on to the rotationally-induccd 
quadruple moment. Rotating configurations of con- 
stant density stars can be described by Maclaurin 
spheroids [89]. The quadrupole moment in Newtonian 
theory is given by [90] 

p-K pr.ie) 

Q(rot),N ^ 27r / / p{r, ey"^ P2{cos 6) sin 9 dr d0 . 
Jo Jo 

(64) 

The surface of the star r ~ r* (9) is in turn given by 



M0) 



sin^ 9 



cos 



62 



'1/2 



(65) 



where b and c are the semi-major and semi-minor axes, 
respectively. By substituting Eq. (65) and p = pc into 
Eq. (64), one obtains 

Q("*''~ = -^Pc&M6'-c')- (66) 
1 o 

From the equation of hydrostatic equihbrium, 

^ = --Vp-VvI/, (67) 
at p 



where boldfaced quantities refer to three-dimensional Eu- 
clidean vectors, with v = ft x r and the NS angular 
velocity vector, one obtains [89] 



2Trpc 



Vr^(3-2e2) . _i 3(1 -e2) 
sm e ; 



1 ^ 1/2 



/Stt 
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Pce + 0{e^), 



(68) 



where 51, = |r2| and e is the eccentricity defined by 



&2 



(69) 



Using Eqs. (68) and (69), we can eliminate c from 
Eq. (66) and substitute 6 = 7^* -I- 0{nl) to obtain 

Q^"*)'" = -^7^^o2 + o(^74). (70) 

Keeping only the leading term, one obtains [90] 

Q-^f^^ (71) 

The dimensionless rotational Love number A'-'^°*^'" can 
be calculated as 



(72) 



which agrees with A^*''^^'" given in Eq. (63). This then 
verifies Eq. (59). 

From Eqs. (62), (63) and (71), one obtains the I-Love- 
Q relations in the Newtonian limit for the incompressible 
EoS: 



(n=0) 
IX 



with 



^(rot) ,N 


2/5 


^(rot) ,N 


1/5 


(n=0) _ 


27/5 


I\ 




(n=0) 


128 


TQ 


3125 


(n=0) _ 


25 


QA 


214/5 



0.528, 
i 0.0410, 
a 3.59. 



(73) 

(74) 
(75) 
(76) 



Of course, universality would be established if the con- 
stants C'a are independent of the EoS, with A any pair 
in the I-Love-Q trio. We will compute the same relations 
for the n = 1 polytrope next, and thus, we will verify the 
degree of universality quantitatively. 



3. Polytropic I-Love-Q relations: n = 1 

Let us now concentrate on the n = 1 polytrope and 
look first at the moment of inertia. Equation (12) in 
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the Newtonian limit gives the equation of hydrostatic 
equihbrium: 



dp 
dR 



pM 
1^ 



(77) 



From Eqs. (10), (77) and p ~ Kp^, one can solve this 
equation to obtain 



1 AL 1 



ttR 



(78) 



the n = and n = 1 polytropic EoSs in GR approach 
the Newtonian ones as the compactness decreases. No- 
tice, however, that we have here analytically shown that 
the I-Love-Q relations are very similar for the n = and 
n = 1 polytropic EoSs only. This does not mean that the 
dependence of the I-Lovc-Q relations on the EoSs is weak 
in the Newtonian limit for all EoSs. Indeed, Fig. 9 shows 
that these relations for some EoSs, such as APR and SLy, 
do not approach the Newtonian limit of the n = and 1 
polytrope. 



One can then calculate (and I^) by substituting the 
above equation in Eq. (25) to find 



2(7r2-6) 1 



37r2 



C2 



(79) 



Let us now consider the tidal apsidal constant. The 
solution to Eq. (60) can be written in terms of Bessel 
functions, as oc (i?/7^*) J5/2(7ri?/7^*) [8, 9]. With 
this, we find that the apsidal constant is 



(tid)^^ 



15 
2^ 



(80) 



This constant agrees with the numerical results of [86 
The dimcnsionless tidal Love number, A^'"^^'", is then 



^(tid), 



15 -TT-^ 1 
3^2 C5 



(81) 



Finally, let us look at the rotationally-induccd 
quadrupole moment. From Eq. (59), one easily finds that 



37r2(15-7r2) 1 
4(7r2-6)2 C 



(82) 



We now have all the necessary ingredients to compute 
the LLove-Q relations in the Newtonian limit for an n = 
1 polytrope. From Eqs. (79), (81) and (82), one finds 



(n=l) 
I\ 



with 



^(n=l) 



^(tid),N 
^(tid),N 



2/5 
1/5 



lA 33/57r6/5(15-7r2)2/5 



0.527, 



C 



(„=1) _ 32(^2 _ Q)5 



IQ 



277r6(7r2 - 15)2 

36/5^12/5(15 _ ^2)4/5 

4(7r2 - 6)2 



0.0406, 



3.60. 



(83) 

(84) 
(85) 
(86) 



Observe that the numbers shown in Eqs. (84)-(86) are 
almost identical to those in Eqs. (74)-(76). 

The I-Love-Q relations for the n = and n = 1 poly- 
tropic EoS in the Newtonian limit are shown as black 
thin solid lines in Fig. 9. Notice that the relations for 



4- Analytical 



The universality of the LLovc-Q relations rests on two 
ingredients. The first ingredient is that the functional 
form of the relations must be the same for different EoSs. 
For the n = and n = 1 polytropes in the Newtonian 
limit, this is verified by comparing Eqs. (73) to (83), 

and noting that regardless of the EoS, oc [A**'^)-"] 

oc [Q~]^ and Q~ oc [A^^^^)'"] This fact is perhaps 
expected, since all multipole moments must be propor- 
tional to the product of a dimcnsionless constant and the 
compactness to some power. The power is determined 
by the Newtonian dimensional structure of the particu- 
lar multipole moment, e.g. / oc 7?,2 and thus / oc C~^. 
This power will be the same regardless of the EoS, and 
thus the power exponent in the LLove-Q relations will 
also be EoS independent. 

The second ingredient, and perhaps the most difficult 
to understand, is the requirement that the constants of 
proportional (ie. the Ca's) be the same regardless of the 
EoS. For the n = and n = 1 polytropes in the Newto- 
nian limit, this is again verified by noting that the coef- 
ficients in Eqs. (74)-(76) are almost identical to those in 
Eqs. (84)-(86); their ratios are 



^(»-0) ^ 22/533/5^6/5(15 _ ^2)2/5 ^ ^ 

c' 



002 , (87) 



IX 

_ 1087r6(7r2 - 15)2 



3125(7r2 - 6)5 

IQ 



1.008. 



25(7r2 - 6)2 



^(n=0) 

^QA ^ 

(j{n=D 24/536/57^12/5(15 _ 7r2)4/5 
QA 



0.997. (89) 



In principle, there is no reason to expect that these 
coefficients should be equal to each other regardless of 
the EoS. Rather, one expects them to depend on the NS 
internal structure. One possible explanation is to argue 
that these coefficients depend on integrals of the energy 
density that are more heavily weighted toward the NS's 
outer layers, i.e. the structure of the NS in its outer layers 
is what is mostly determining these coefficients. But it is 
precisely in the outer layers that nuclear physics uncer- 
tainties are lowest. Therefore, the EoSs in this regime are 
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more similar to each other than in the core, thus leading 
to some degree of universality. 

We have found some evidence to support this interpre- 
tation, shown in Fig. 9. Focus on the I-Love-Q relations 
for the polytropic EoSs with n = 2, 2.5 and 3, which 
greatly modify the internal structure in the NS's outer 
layers, far from the core. In fact, these poly tropes lead 
to essentially no energy density near the NS surface, with 
most of it concentrated near the core. We see that, in- 
deed, when we choose an EoS that affects the NS struc- 
ture far from the core, we significantly lose universality in 
the I-Love-Q relations, as one can see from the n = 2, 2.5 
and 3 curves in the bottom part of each panel in Fig. 9. 




FIG. 11. (Color Online) {p/ pc){R/'R.,f (top) and 
p{dp/ dp) / {pj./pc) (bottom) as functions of R/TZ, for C = 0.17 
with various EoSs. The former corresponds to the integrand 
(modulo normalization constants) of J" [Eq. (24)] and Q^"*' ''^ 
[Eq. (64)] in the Newtonian limit, while the latter corresponds 
to the EoS-dependent coefficient in Eq. (60), which gives 
A'""^''". Sudden chan ges in the bottom panel correspond to 
nuclear phase transitions. Observe that the dominant contri- 
bution of (p/pc)(J?/7?.*)* and p{dp/dp)/{pi/pc) for realistic 
EoSs come from the NS outer layer. 

Further evidence can be found by investigating a few of 
the terms that control the behavior of the moment of in- 
ertia, the quadrupole moment and the tidal Love number 
as a function of C. First, in the Newtonian limit, both 
[Eq. (24)] and [Eq. (64)] can be written in inte- 

gral form, where the radial dependence of the integrand is 
proportional to p{R)R^. The top panel of Fig. 11 shows 
{p/pc){R/n^Y as a function of i?/7^, for C = 0.17 with 
various EoSs. / and Q'"'"') are proportional to the area 
under the curves in this panel. Observe that the curves 
are similar, even for different realistic EoSs, and the dom- 
inant contribution comes from the NS outer layer, some- 
where between R/TZ^, w 0.7 — 0.9. This may explain the 
similar behavior of the I-C and Q-C curves with different 
realistic EoSs in the high-compactness regime (see Figs. 6 
and 7), as well as the universal I-Q behavior. Observe 
also that as one increase the polytropic index n, the NS 



becomes more centrally-condensed, and / and Q arc not 
dominated just by the NS outer layers. 

Similarly, we can study the behavior of the structure- 
dependent term that determines the tidal Love num- 
ber in the Newtonian limit. Equation (60) shows that 
there is only one such term and it is proportional to 
pdp/dp — p{dp/dR){dp/dR)~^ . The bottom panel of 
Fig. 11 plots p{dp/dp)/{p1lpc) as a function of R/TZ>, 
for C — 0.17 with various EoSs. Sudden changes in the 
slope corresponds to nuclear phase transitions. Similar 
to the top panel, the behavior of pdp/dp is similar among 
realistic EoSs and the dominant contribution comes from 
the NS outer layers. This partially explains the similar 
behavior observed in the Love-C curves with different re- 
alistic EoSs in the high-compactness regime (see Fig. 8), 
as well as the universal Q-Love and I-Love relations. 

Another possible explanation for the universality of the 
I-Lovc-Q relations involves the behavior of this trio as 
one approaches the BH limit. The no-hair theorems [43- 
48] of GR state that the exterior multipolar structure of 
an isolated, stationary, axisymmetric BH solution in GR 
is completely determined by its mass and its spin angu- 
lar momentum. Therefore, the quadrupole moment, for 
example, is completely determined by the spin angular 
momentum through Eq. (34) with A = [49, 50]. For 
NSs, such a result docs not exist, but one might still ex- 
pect the TQ relation to become less structure dependent 
as one approaches the BH limit (C — > 0.5). Indeed, Fig. 9 
shows that the loss of universality (the relative fractional 
difference shown in the bottom part of each panel) de- 
creases logarithmically as the mass increases from O.3M0 
to 1.4M0. 

Of course, one can never increase the compactness 
enough by a finite amount to turn a NS into a BH, ie. the 
NS sequence of varying compactness does not terminate 
in a BH for finite central density. Still, it is interesting to 
see that as C increases, the I-Love-Q relations are indeed 
approaching the BH limit, as explicitly shown in Fig. 9. 
Moreover, universality in the I-Q curve suggests a univer- 
sal relation between the NS spin and the NS quadrupole 
moment that is almost independent of the internal struc- 
ture. Such a relation is similar to the no-hair relations 
for BHs [49, 50]. 

Before proceeding, let us point out that this effacing 
of internal structure is not the same as what is discussed 
in the effacement principle [51] in GR. The latter states 
that the equations of motion of compact objects of any 
size and structure depend only on integral parameters, 
like the mass and spin, and it is independent of its actual 
shape and internal structure. Of course, this principle 
holds in GR but only for BHs because of the no-hair the- 
orems. The effacement principle is violated for NSs, but 
the violation is small, with corrections to the accelera- 
tion entering at 5PN order for a binary of non-spinning 
compact objects. Since the effacement principle deals 
with the motion of a body only, and not on the multipo- 
lar structure of its exterior gravitational field, the efface- 
ment we find here is not a consequence of the standard 
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efFacement principle. 



VIII. APPLICATIONS 

The I-Lovc-Q relations have 3 immediate applications 
to observational astrophysics, GWs and fundamental 
physics. Let us look at each application in turn. 



A. Observational Astrophysics 

On an astrophysical front, a measurement of any sin- 
gle member of the I-Love-Q trio automatically provides 
information about the other two, even when measuring 
the other two directly might not be possible with current 
observations. For example, one might be able to measure 
/ within 10% accuracy by measuring the orbits of binary 
pulsars sufficiently accurately, so as to extract the spin- 
orbit coupling effect in the advance rate of the periastron 
of the double binary pulsar J0737-3039 [6, 7]. If such 
measurement is accomplished, one can then automati- 
cally obtain the quadrupole moment and tidal Love num- 
ber of the primary pulsar by using the I-Love-Q relations. 
Similarly, if an equal-mass NS binary within 300Mpc is 
about to coalesce, one might be able to determine the 
^(tid) p£ ^YiQ constituents with second-generation, ground- 
based GW interferometers [17-19, 25, 26]. Then, from 
the I-Love-Q relations, one may obtain the moment of in- 
ertia and quadrupole moment of the binary constituents, 
which again, would be very difhcult to measure with 
GWs. 

Of course, these application assumes that the univer- 
sality of the I-Love-Q relations holds, which rests on the 
assumptions of uniform and slow-rotation, small tidal 
perturbations and that GR is the correct theory. There- 
fore, this technique cannot be applied to, for example, 
newly-born NSs that are differentially and rapidly ro- 
tating. NSs that source GWs in the sensitivity band of 
ground based detectors, however, are expected to be old, 
and thus uniformly rotating with large spin periods (they 
should have spun down by the time they are visible by 
GW detectors [60]), so that the slow-rotation approxi- 
mation is well-justified. The primary NS in the double 
binary pulsar has a period of 22 ms, which implies a 



X ^ 0.018, small enough that the slow-rotation approxi- 
mation is again well-justified. 

Millisecond binary pulsars with short-periods, ie. peri- 
ods below 1 ms, would be spinning too fast for the above 
relations to be directly applicable. However, we expect 
I-Love-Q type universality with respect to the EoS to 
still hold in this case, except that now the coefficients 
in Table I will also depend somewhat on the spin angu- 
lar frequency. One can correct the universal I-Love-Q 
relations for non-negligible spin by considering rapidly 
rotating NSs [61-64], but wc leave this to future work. 
B. Gravitational Wave Astrophysics 

Another application of the I-Love-Q relations is to GW 
astrophysics, as a means to break the degeneracy between 
individual spins and the quadrupole moments of NSs in 
the GWs emitted during binary NS inspirals. Let us 
first discuss gravitational waveforms of spinning, tidally- 
deformcd NS binaries, and then, carry out a back-of-the- 
envelope parameter estimation study using Fisher theory. 
The latter will allow us to determine the degree to which 
degeneracies arc broken through the I-Love-Q relations 
and the projected accuracy to which individual NS spins 
could be measured given a GW detection. 



1. Waveforms 

The sky-averaged gravitational waveform (in the 
Fourier domain) generated by a compact NS binary in 
a quasi-circular orbit with masses mi and m2 and at dis- 
tance Dl is given by [91] h{f) = A(/) cxp[i*(/)], with"* 



A{f) 
*(/) 



1 



7\/(5/6 



V30^2/3 Dl 

*tp(/) + *Q(/) 



%(/). 



(90) 
(91) 



Here, / is the GW frequency, A4 = mrf^^ is the chirp 
mass, 77 = mim2/iTi^ is the symmetric mass ratio and 
TO = 7711 + 7772 is the total mass. The quantity Vl'tp(/) 
is the gravitational waveform phase in the test-particle 
limit while and represent terms that deviate from 
this limit, where the former corresponds to a quadrupole 
moment deformation, while the latter depends on the 
tidal Love number. 

The test-particle term, to 3.5 FN order, is given by [92- 
94] 



*tp(/) = 2nft., - 0c - 5 + ^^i-Mfr^/'l l+^^ + fr^y- (167r - 4/3)x^/2 



/ 15293365 27145 3085 , „ \ 2 / 38645 65 \, ^, ^ 

1 7? H ir - 10a hz: + TT nri - 7 (1 + 3 log v)x^''' 

\ 508032 504 ' 72 ' / I 756 9 ' " ^ ^ ' 



* needs to be multiplied by •^3/2 when calculating the Fisher 

matrix for LISA and DECIGO/BBO. 
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where X = ~ (Trm/)^/'^ and {tc, (l>c) correspond to the 
time and phase at coalescence, respectively, with 7e the 
Euler constant. The spin parameters /3 and cr [55, 93, 95], 
7 [93] and a [94] are given by^ 



113 19 , , ^ 



(93) 



719 , 
1^^ 



L ■ Xa 



719 1 , , , 



233 

Is" 



^m{Xs ' Xa) 



on 



_ '233 T_ 
^ ' "96" 24^^ 



/ 732985 
V 2268 ~ 
/ 732985 



Xs 

24260 



V 2268 
22707r ' 



a 



1 



81 
140 
~9~ 
227 
156 



ri^ S,n(L- Xa 

f]] (l • xJ + S, 



L Xc 



L ■ Xa^ 

(94) 

(95) 

(96) 



where L is the unit orbital angular momentum, 8m = 
(mi — m2)/m is the dimcnsionless mass difference, Xs = 
iXi + X2)/2 and Xa = ixi - X2)/2 with Xi = S^lm\ 
denoting the dimcnsionless spin vector of the i-th body. 
Notice that we are here referring to the individual NS 
spin vectors by Si 

The quadrupole moment contribution correction to the 
test-particle limit in the GW phase enters at 2PN order 
and it is given by [80, 95] 
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Qs is strongly correlated with cr, which enters at the same 
PN order as Qs- 



The leading-order contribution of to the GW 

phase enters at 5PN order through [17] 



= -l|^24 [(1 + 7r; - 31r,^)A, 



where 



T(tid) T(tid) , ^ 

X _ ^1 + ^2 \ — '^l 
= -Z , Aa = 



(tid) T(tid) 
" ^2 



(99) 



(100) 



Higher PN contributions to "^xif) can be found in [23, 
24, 26]. When carrying out parameter estimation studies, 
as explained below, we will use '^\{f) as given in [26], 
which includes up to 2.5PN order corrections relative to 

yJ;5PN6 



2. Parameter Estimation 

For stationary and Gaussian detector noise, the mea- 
surement accuracy of parameters 0" can be estimated as 



(r-^)° 

N 



(101) 



where A'^ is the number of effective interferometers and 

(102) 



^'L SM) 



is the Fisher matrix, where the partial derivatives are 
with respect to the parameters 9°'. The noise spectral 
density 5'„(/) is given in Refs. [96-98] for Adv. LIGO, 
ET and DECIGO/BBO, respectively. We take the lower 
cutoff frequencies /min to be lOHz for Adv. LIGO, IHz 
for ET and the frequency lyr before coalescence for DE- 
CIGO/BBO. For Adv. LIGO and ET, we take the higher 
cutoff frequency to be that of the innermost stable cir- 
cular orbit (ISGO), /n,ax = /isco = 1/(63/2^), while 
we set /max = lOOHz for DECIGO/BBO. N is the num- 
ber of effective interferometers, which we take to be 5 for 
2nd-generation ground-based detectors (corresponding to 



^ cr includes the quadrupolc-monopole interaction in the tcst- 
particle hmit. 



® Tidal effects on the gravitational waveform phase have been cal- 
culated to 1.5PN order relative to the leading 5PN contribution, 
in addition to tail effects at 2.5PN order. Ref. [26] estimated 
that currently unknown terms should be subdominant, at least 
for an equal-mass binary. 
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2 Adv. LIGO, Adv. VIRGO, KAGRA and INLIGO), 2 
for ET (like LISA [99]) and 8 for DECIGO/BBO [98]. 

We focus here on GWs emitted during the quasi- 
circular inspiral of NSs with aligned spins, since this is a 
realistic astrophysical scenario [100]. Given that the NS 
masses are expected to be approximately the same, we 
will not include Qa and \a in the parameter vector, as 
this must be close to zero. We will consider the case of 
slightly unequal NS masses. We choose two parameter- 
izations of the waveform. Parameterization A uses the 
parameter set [26] 

{e\} = (ln7W,ln77,/3,7^i,i„0„A,) (103) 

with the priors |r/| < 0.25 and |/3| < 0.2. We do not in- 
clude (J in this set because the NS spins at the time of 
coalescence are expected to be small [60] ' . Paramctriza- 
tion B uses the parameter set 

{e's} = i\nM,6,n,Xs,Xa,DL,tc,cl)c,QsCXs),~Xs), (104) 

with the priors \Sm\ < 1/3, \xs\ < 0.1 and |xa| < 0.1. 
Moreover, we use the Q-Love relation to express Qs in 
terms of As, and thus, partially break the degeneracy 
between Qs and Xs- 

Figure 2 shows the measurement accuracies of spin pa- 
rameters using second-generation, ground-based detec- 
tors. We assume that that the detected GW was emitted 
by a source at Dl = lOOMpc with SNR - 30. We con- 
sider 3 different systems: (i) (mi, 7712) = (1.45, 1.35)M0, 
XI = X2, (ii) (TOi,m2) = (1.45,1.35)Mq, xi = 2x2 
and (iii) {mi, 1712) = (1.4, 1.35)M0, xi = X2- Observe 
that A/3 and Axa in Fig. 2 are almost identical for all 

3 systems; this is because the measurement accuracy for 
these parameters is dominated by their priors. Observe, 
however, that the averaged spin Xs can be measured to 
0(0.01). Such an accuracy on Xs is inaccessible without 
the Q-Love relation. 

We can understand this enhanced accuracy in the ex- 
traction of Xs as follows. First, notice that, for an 
equal-mass and spin-aligned binary, f] ~ C'(10)xs [see 
e.g. Eq. (93)]. Given that the measurement accuracy of 
/? is Ap = C(O.l), this implies a measurement accuracy 
of Xs of Axs ~ 0.01. The measurement accuracy of Xs 
for system (ii) increases as xi ~^ 0.1, as shown in Fig. 2. 
This is because dh/d6m ~ w dh/dxa as xi ~ 0.1 and 
the priors lead to <5m, Xa and other parameters being 
effectively uncorrelatcd. In such a case, however, the as- 
sumptions that underlie the Fisher approximation may 
be violated [101], and hence, one requires a Bayesian 
analysis [96] to confirm these results. Finally, we have 
checked that the Q-Love relation does not improve the 
measurement accuracy of As . 

Up until now, we have considered equal-mass NS bina- 
ries, but realistic systems may not have identical masses. 
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FIG. 12. (Color Online) Separatrix between systems with 
'^x^ > 1 (above the curves) and ^^^^^ < 1 (below the curves) 
as a function of Am = jmi — m2| and fh = (mi +m2)/2, and 
for different realistic EoSs. For systems below the curves, we 
can safely neglect Xa in parameter estimation, provided the 
SNR « 0(10). is smaller than for |x| < 0.1. 



If the NS masses are not equal, one must then take into 
account the parameters Qa and Xa , whose inclusion could 
in principle degrade the accuracy to which other parame- 
ters are extracted. Let us then study the range of masses 
for which neglecting Qa and Xa is a good approximation. 
A rough estimate of this range can be obtained by in- 
vestigating the systems for which the accumulated GW 
phase induced by terms proportional to Qa and Xa is less 
than one radian. Let us then define the latter by 
and 5*3^ respectively, where 



75 1 
64^ 
9 1 



4x? 



2 ,,2 



xi Qa 



-1/2 



2\\ X 2,5/2 



(105) 



16 ?7 

to leading PN order. 

Figure 12 shows the range of masses for which "^f^^ = 1 
for various realistic EoSs. Systems above these lines 
would lead to ^'^ > 1, while those below this line lead to 
< 1. In particular, for systems that satisfy the lat- 
ter inequality, we can in principle neglect Xa if the SNR 
is 0(10). This figure implies that for second-generation, 
ground-based detectors, the parameter estimation study 
presented above is probably valid, even for unequal-mass 
systems provided, for example, that ttt. = 1.4Mq and 
A777 < O(O.1)M0. where 7fis = (7711 -|- 7772)/2 is the av- 
eraged mass of the binary and A777 = [7771 — 7712] is the 
mass difference. This same conclusion also applies to ne- 
glecting Qa, provided [x] < 0.1. 



C. Fundamental Physics 



Damour et al. [26] estimated that at the time of coalescence, 
|/3| < 0.2 and \a\ < 10"''. 



The independent measurement of any 2 members of 
the I-Love-Q trio would allow us to perform model- 
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TABLE II. NS parameters for A("<^' = 400.0 and M, = 
1.3382Mq. /spin corresponds to the NS spin frequency with 
X = 0.1. 

independent and EoS-independent tests of GR. For ex- 
ample, if one can measure / and A*-*"^-* independently, one 
can plot a point in the I-Love plane with an error box. 
If the I-Love relation in GR crosses the error box, then 
GR is consistent with the observations. Otherwise, one 
would have found mo del- independent evidence for some 
type of departure from GR. Moreover, one can constrain 
non-GR theories by requiring that the I-Love relation in 
that theory crosses the error box. 

The accuracy of such a test depends, of course, in how 
accurately two elements in the I-Love-Q set can be mea- 
sured. One way to measure / would be to look for a 
spin-orbit correction to the rate of advance of the peri- 
astron of a binary system. Future double binary pulsar 
observations may measure / with an accuracy of roughly 
10% [6, 7]. Probably, the best way to measure A^'"^' and 
Q would be to use GW observations. 

In what follows, we first discuss the possibility of mea- 
suring Q and A*^*"^^ simultaneously, given GW observa- 
tions. Then, we discuss how well GR tests can be carried 
out by combining GW observations with binary pulsar 
observations. For concreteness, we apply all of this to 
a specific modified gravity theory (dynamical CS grav- 
ity [57]). 

1. Redundancy Tests with GW Observations Only 

Let us estimate how accurately future ground-based 
detectors may determine A^*"^-* and Q simultaneously. 
The measurability of A*^*'^-* has been discussed extensively 
in [17-26], while the effect of the quadrupole moment on 
compact binary GWs has been estimated in [60, 80, 102], 
but so far no study has been performed to study their si- 
multaneous extraction. Let us consider an equal-mass 
NS binary with A^*"^^ = 400, which corresponds to NSs 
with = 1.4M0 with the APR EoS (other parameters 
are shown in Table II). 

Figure 13 shows the measurement accuracy of Qs and 
As with second-generation, ground-based detectors, ET 
and DECIGO/BBO. We assume an equal-mass, spin- 
aligned NS/NS binary with (1.4, 1.4)M0, xi < 0.1 and 
X2 = at = lOOMpc, and we also assume that the 
APR EoS is the correct one. The measurement accuracy 











— LIGO 






- • ET 






■ - DECIGO/BBO 

















FIG. 13. (Color Online) Measurement accuracies for Qs 
(top) and Xa (bottom) given a GW detection, emitted by 
a (1.4,1.4) Mq, spin-aligned NS/NS binary system with 
Adv. LIGO, ET and DECIGO/BBO. We assume = 400, 
X2 = 0, Di = lOOMpc and the APR EoS, with priors 
\S\ < 1/3, |xi| < 0.1 and |x2| < 0.1. Although it may be 
difficult to measure Q,, it should be possible to measure As 
with ground-based detectors. 



of AQs with second-generation, ground-based detectors 
is AlnQs ~ 30, which would increase to AlnQs « 2 
using future detectors, such as ET or DECIGO/BBO. 
These results imply that it may be difficult to measure 
AQs due to its strong degeneracies with spin parame- 
ters. Notice, however, that even though one may not 
be able to detect Qs, one can still place an upper and 
lower bound on Qs- Such a bound would be sufficient to 
perform model-independent GR tests. The measurement 
accuracy of As with second-generation, ground-based de- 
tectors is AAs ~ 0.8, which would increase by roughly an 
order of magnitude using ET. Although the error bars 
arc large with current detectors, it may be possible to 
measure As with future GW observations. 

Given the above measurement errors, we can now simu- 
late a GR test. Figure 14 presents the Q-Love relation for 
realistic EoS, together with a fiducial GW measurement 
of the pair (Q, A'-'"^-') and its estimated errors. Notice 
that the error in Q is larger than the value about which 
the error is centered. This implies that a GW measure- 
ment would not be able to measure Q, but it would be 
able to say the region of allowed Q that is consistent 
with the GW detection. Therefore, such a GW detection 
would automatically constitute a model-independent test 
of GR; for GR to be consistent with these measurements, 
the GR Q-Love curve must cross the GW error box in 
(Q,A("d)). 

The above test is quite robust. First, although we 
employ a uniform and slow-rotation approximation, the 
NS /NS binaries that ground-based detectors will observe 
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FIG. 14. (Color Online) The Q-Love relations for realistic 
EoSs with proposed measurement errors. We assume that 
we determine Q and A''"^-* simultaneously by detecting GWs 
from an equal- mass, spin-aligned NS/NS binary with xi ~ 
0.1, X2 = and A^*'"*' = 400 at Dl = lOOMpc with ET. 
The fiducial values of (A^'''^^ Q) are shown as a big black 
cross. Although the measurement error of Q is greater than 
its fiducial value, one can still perform model-independent 
and EoS-independent tests of gravity by requiring that the 
Q-Love curve must pass through the error box. 



will have spun down by the time they enter the detector's 
sensitivity band, and thus, the slow-rotation approxima- 
tion should be excellent. Second, the error box of Fig. 14 
depends on how accurately (Xs,Qs) can be measured, 
which in turn depends on whether (Xa,Qa) need to be 
included in the parameter set. This would be the case if 
the binary system detected were not an equal-mass one. 
As shown in Fig. 12, however, there is a wide range of 
mass ratios for which these parameters can be neglected, 
even outside of the equal-mass point; thus, the discussion 
presented above should be robust. 



2. Joint Tests with GW and Electromagnetic Observations 

Since Q is a quantity that is difficult to measure with 
GW observations, let us consider model-independent and 
EoS-independent tests of GR with the TLove relation 
that uses a combination of GW and double binary pul- 
sar observations. Let us then assume that / has been 
measured to 10% by future double binary pulsar (J0737- 
3039) observations [6, 7], and that A''"^-' has been mea- 
sured to 60% with future GW observations. The latter 
assumes an ET detection of an equal-mass, non-spinning 
NS/NS binary at 3Gpc, with the individual NS masses 
exactly equal to that of the primary pulsar in J0737-3039, 



= 1.3382M0, assuming the Shen EoS®. All of this is 
shown in Fig. 3, together with the fiducial measurement 
of (/,A("'i)) as a big black cross. As shown in that fig- 
ure, one can constrain modified theories of gravity, such 
as dynamical CS gravity, by requiring that the TLove 
curve crosses the error box. 

The test described above has one major problem: the 
NS mass rripuisar of the primary pulsar in J0737-3039 
and the individual NS masses in the binary system that 
generated the detected GW will all in principle be dif- 
ferent from each other. As explained in Sec. VIII B 2, 
the accuracy to which As can be measured assumed that 
Qa and Aa could be neglected, which holds for a cer- 
tain range of mass differences Am, shown in Fig. 12. 
In general, the typical maximum mass difference would 
need to be Am — O{0.1)Mq, assuming observations with 
SNR « 10. Given current event rate estimates, one ex- 
pects to detect NS/NS binaries with such similar masses, 
and thus, this is not in principle a problem. The test 
above, however, also requires that TOpuisar ~ 'fn, since af- 
ter all the TLove-Q relations assume one is investigating 
NSs with the same mass. 

Let us then estimate how much the I-Love-Q rela- 
tions would change if rripuisar 7^ rrrGw- The top panel of 
Fig. 15 shows the I-Love relation with rripuisar/™Gw = 
0.9, 1.0 and 1.1 for realistic EoSs, while the bottom 
panel shows the relative fractional difference between 
the I-Love curves for different EoSs and the APR EoS 
as a reference. The relative fractional difference when 
TOpuisar/rrrcw = 0.9 (not shown in this figure) is similar 
to that of rripuisar /rrr-Gw = 1.1. One sees that the depen- 
dence on the EoS becomes stronger as the mass difference 
riTpuisar and rrrcw increases. However, this dependence is 
still weak if the mass difference is sufficiently small (or 
order O.IMq), and most importantly, the loss of univer- 
sality (the difference between curves with different EoS) 
is much, much smaller than the observational error in 
measuring either the moment of inertia or the Love num- 
ber. Therefore, one can perform the GR test described 
above, even when TOpuisar 7^ TOgw 

Of course, the test described here assumes that the uni- 
form and slow-rotation approximation used to derive the 
I-Love-Q relation holds for binary pulsars. This is indeed 
the case, provided the period is sufficiently long, such 
that each binary component is spinning slowly. However, 
the approximation might break down for (currently un- 
observed) sub-millisecond pulsars, ie. those with periods 
shorter than 1 ms. For such systems, the I-Love-Q rela- 
tions will also now depend on the spin frequency. A cur- 
sory analysis, however, suggests that the spin-frequency 



* Adv. LIGO is expected to detect NS/NS binaries out to Dj^ Ri 
SOOMpc with the detection rate of C'(10)/yr [13]. Therefore, if 
we consider ET detecting GW signals from a NS/NS binary at 
3Gpc, the expected detection rate would be C'(10)/yr X 10'^ ~ 
C'(10*)/yr. With this detection rate, we may detect an equal- 
mass NS/NS binary, with the individual masses very close to 
that of the primary pulsar of J0737-3039, M. = 1.3382M0 
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FIG. 15. (Color Online) (Top) 1-Love relation of NSs with 
mass ratios of mpuisar/oT-GW ~ 1 (solid), 1.1 (dashed) and 0.9 
(dotted-dashed) for realistic EoSs. (Bottom) Relative frac- 
tional difference with the APR curve as the reference. Ob- 
serve that the loss of universality when rUp^isai / rriGw = 0.9 is 
similar to that when mpuiaar/wGw ~ 1.1- Observe also that 
the loss of universality is small relative to the observational 
error in measuring the moment of inertia or the tidal Love 
number. 



effect breaks universality at the 10% level [61-64]. There- 
fore, the difference in I-Love-Q relations for different 
EoSs will be rather small, and in particular smaller than 
the errors in the first binary pulsar measurements of the 
moment of inertia. 

Binary pulsar observations may measure / within 10% 
accuracy, but this amazing measurement will be difficult 
to accomplish in the near future [7] . This is because the 
effect of the moment of inertia (or equivalently, the spin- 
orbit coupling) in the motion of the binary is of 0(1O~^) 
relative to the leading-order contribution. This means 
that one needs to measure at least 3 post-Keplerian pa- 
rameters to this accuracy, in order to determine the two 
masses and the moment of inertia. Of course, such a 
measurement is a big challenge, although not out of the 
question, as binary pulsar observations improve within 
the next ten years. 

Alternatively, one could use the NS compactness C 
instead of / to perform model-independent GR tests. 
Currently, C has been measured to C(10)% accuracy 
with low-mass X-ray binary observations [2-5]. Figure 16 
shows the A^*"^-'-C relation for realistic EoSs and for the 
n = 1 polytrope. We also show in this figure a fidu- 
cial measurement of (C, A^*"^)) (big black cross), as well 
as projected measurement accuracies (dashed lines). For 
the latter, we assume AC = 0.05 from electromagnetic 
observations of a NS with A/* = IAAIq and a roughly 
70% measurement of A^'"^^ from GW detection (equal- 
mass, non-spinning NS/NS binary with the individual 
NS mass of Ah = 1.4A/0 at I^l = lOOMpc) with second- 



generation, ground-based detectors^. Although the de- 
pendence on the EoS is relatively large compared to the 
universality of the I-Love-Q relations, the measurement 
errors are larger than the uncertainties due to the EoS. 
This shows that one might be able to use the Love-C 
relation to perform model-independent tests of gravity. 




c 

FIG. 16. (Color Online) Love-C relations with realistic EoSs 
and the n = 1 polytrope (solid lines), a fiducial measurement 
of (A'^'^^C) (black cross) and projected measurement uncer- 
tainties (dashed black lines). We assume AC = 0.05 for a NS 
with mass 1.4M0, and a measurement of A with a roughly 
70% error. Observe that the Love-C relation loses some of 
the universality shown in the I-Love-Q relations. However, 
the error introduced due to EoS dependence is much smaller 
than the measurement error in the compactness or the tidal 
Love number. 



3. Example: Dynamical CS Gravity 

We now apply the results obtained above to sec how 
testing a specific theory of gravity would go about. As 
an example, we choose dynamical CS gravity [57, 103], 
which is well-motivated from the Standard Model, su- 
perstring theory [104, 105], loop quantum gravity [106- 
108] and infiation [109]. Dynamical CS gravity is a 
parity- violating, quadratic-curvature theory, where the 
Einstein-Hilbcrt action is modified through the Pon- 
tryagin density (the contraction of the Ricmann ten- 
sor and its dual), coupled to a dynamical scalar field. 
This theory has a characteristic length 5}^^, which has 
been constrained by Solar System experiments, using 
Gravity Probe B [110] and LAGEOS [111], to S}'"^ < 
O(10^)km [56]. Dynamical CS gravity should be treated 
as an effective theory, and thus, one should work to 



^ The measurement error of AA^'"*) 0.7 is slightly better than 
that shown in Fig. 14. This is because we assumed parameter 
set A [Eq. (103)] instead of B [Eq (104)]. 
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leading-order in a small coupling expansion, ie. to lead- 
ing order in the dimensionless coupling constant C = 

man [67]. 

NSs in dynamical CS gravity have been studied before. 
Reference [67] found that it would be difficult to mean- 
ingfully constrain this theory with binary pulsar observa- 
tions in the standard fashion. This is because the largest 
CS correction appears in the rate of change of periastron 
advance at IPN order, and thus it is suppressed by the ra- 
tio of the binary's mass to its separation (for J0737-3039, 
this is of 0(10"^)). The CS correction to the NS moment 
of inertia was calculated in [56, 73], while the CS correc- 
tion to the NS quadrupole moment was obtained in [67] . 
In the small coupling approximation, the CS corrections 
to / and Q scale linearly with C^. The leading-order CS 
correction to tidal effects enters through the gravitomag- 
netic tidal tensor (because of the parity) [59] , and hence 
^(tid) affected at leading-order. 

Figure 3 shows the I-Love relation in dynamical CS 
gravity with a fixed value of the coupling constant ^ = 
1.1 X lO^M^. Observe that the dependence on the EoS 
is stronger than that of the OR I-Love relation. We 
believe that this is because a compact object in dy- 
namical CS gravity depends on the scalar-dipole charge, 
which encodes information on the internal structure of 
the body [67]. Given this reasoning, we expect that the 
I-Love-Q relations should be more sensitive to the NSs' 
internal structure in dynamical CS gravity than in GR. 
The bottom panel of Fig. 3 shows that there is indeed a 
loss of universality, but the latter is still preserved to a 
few % level. 

With this in hand, let us estimate the projected bound 
that one could place on dynamical CS gravity using the 
I-Love relation. With ^ = 1.1 x lO^M^, the I-Love curve 
in dynamical CS gravity barely crosses the error box. 
Since the larger ^, the higher the CS curves, such an I- 
Love observation would automatically constrain ^ < l.lx 
lO^M^, which corresponds to C = 0.0581. Converting 
back to dimensional quantities, such a test would impose 
the constraint 

^1/'* < C'(50)km. (106) 

NS observations would then allow us to probe the theory 
within NS length scales, like the NS radius. Notice that 
the above bound is stronger than Solar System [56] and 
table-top [58] ones by more than six orders of magnitude. 
Notice, however, that this bound is slightly weaker than 
the proposed projected bound with GW observations of 
BH/BH binaries [112]. This is because the "radius" of a 
BH is smaller than that of a typical NS, and thus, with 
the former, we can probe shorter length scales. Notice 
also that the bound given above is dominated by the 
measurement error on the NS moment of inertia. This is 
reasonable because the tidal Love number is unaffected 
in dynamical CS gravity to the order of approximation 
considered here. 

The measurement accuracy shown in Fig. 3 is obtained 
by assuming that the Shen EoS is the correct one, which 



gives the weakest bound on the theory among the realistic 
EoSs considered in this paper. This is because, with 
the NS masses fixed to M* = 1.3382M0, the Shen EoS 
gives the largest A^*''^) (see Table II). Figure 3 shows that 
the deviation away from GR becomes larger for smaller 
^(tid)^ This is because the compactness becomes larger 
for smaller A'-'"^-' (or smaller / and Q) which allows us 
to probe stronger gravity. These studies suggest that the 
I-Love-Q relations can be very powerful in testing GR in 
the strong-field regime. 

IX. FUTURE DIRECTIONS 

We have derived relations between the moment of in- 
ertia, the quadrupole moment and the Love numbers, 
I-Love-Q relation, that are essentially independent of the 
EoS for uniformly and slowly-rotating NSs. These re- 
lations open the door to exciting applications in astro- 
physics, GW theory and fundamental physics. We have 
here carried out a preliminary study of a few applica- 
tions, but our paper enables a lot more work. One ex- 
ample is a more detailed study of the measurement ac- 
curacy of binary parameters given a GW observation. 
We here carried out a Fisher analysis, but this is known 
to be inaccurate for signals with low SNR [101], as ini- 
tially expected with second-generation, ground-based de- 
tectors. Such an analysis could be improved on through 
a Bayesian study [96, 101]. Another example is to repeat 
the analysis that uses the I-Love-Q relations to test GR 
to include systems with different mass ratios. This ex- 
tension is particularly important, given that millisecond 
binary pulsars will probably not have exactly the same 
mass as the NSs observed through GWs emitted in the 
late inspiral. Our results suggest that the conclusions 
regarding tests of GR should be robust even when the 
masses differ by 10%, but a more detailed analysis would 
be desirable. 

The analysis presented here has a few caveats that 
should be re-iterated here for completeness, although we 
have discussed this to some extent already in the Intro- 
duction. The framework in which the I-Love-Q relations 
have been found to be essentially EoS independent is one 
that employs a uniform rotation, slow-rotation and small 
tidal deformation approximation. Newly-born NSs are 
expected to be differentially rotating at very short pe- 
riods, where the slow-rotation approximation would not 
be appropriate. Moreover, such NSs are expected to be 
much hotter than those in millisecond binary pulsars and 
those that emit GWs in the band of ground-based detec- 
tors. Temperature could introduce further deviations in 
the universality relations described here. 

Of these limitations, the slow-rotation approximation 
is perhaps the most severe, although such an approxi- 
mation is reasonable for NSs in millisecond pulsars with 
periods comparable or larger than 1 ms. We expect the 
EoS universality found here to persist even when includ- 
ing fast rotation, except that now there will be differ- 



26 



ent universal relations for stars with different spin peri- 
ods. The variation, however, should not exceed 10% [61- 
64] . A possible extension of this work would be to refine 
the universal relations to allow for rapidly spinning NSs. 
This could be achieved by numerically solving for rapidly 
rotating NSs and then extracting its multipole moments, 
as recently discussed in [61-64]. Let us reiterate, how- 
ever, that for almost all NS that have been astrophys- 
ically observed, the spin period is sufficiently long that 
the slow-rotation expansion is an excellent approxima- 
tion. From an academic point of view, however, it is 
also worth studying how differential rotation [113] would 
change the I-Love-Q relations and their universality. 

Another possible extension of our work would be to 
consider more generic NSs with anisotropic pressure [114] 
and large internal magnetic fields. Recent work has 
suggested that in fact the NS interior might be super- 
conducting and super-fluid (see e.g. Refs. [115, 116] and 
references therein). The inclusion of these effects will cer- 
tainly affect the EoS, but it is not clear that this will mod- 
ify the TLove-Q relations presented here. This would be 
particularly so if the I-Love-Q relations are truly only sen- 
sitive to the EoS far from the core, where super-fluidity 
and super-conductivity play a small role. 

Although we have investigated how the I-Love-Q rela- 
tions change in dynamical CS gravity, it would be worth- 
while to study such relations in other modifled theories, 
such as scalar-tensor ones [117] and Einstein-jEther the- 
ory [118]. Given any modified theory, one could investi- 
gate how the I-Love-Q relations change, whether univer- 
sality still holds, and how strong one can constrain other 
theories with future observations. 

An interesting avenue to pursue would be to study 
whether universal relations exist between higher-order, 
multipole moments of the exterior gravitational field of 
isolated NSs. For BHs, the no-hair theorems guaran- 
tee that BH multipole moments can be written entirely 
in terms of the BH mass and spin angular momentum 
(assuming the charge is zero), leading to a unique rela- 
tion that, of course, is independent of internal structure 
(BHs lack any). For NSs, such a relation does not exist, 
since the no-hair theorem does not apply. In this paper, 
however, we have found an interesting relation between 
the quadrupole and the dipole moment of the exterior 
gravitational field of an isolated NS that seems almost 



independent of the NS's internal structure. One might 
then naturally wonder whether similar relations hold for 
higher-order multipole moments, which may lead to a NS 
no-hair conjecture, ie. that NS multipole moments can be 
effectively expressed only in terms of the NS mass M^,, 
the NS angular velocity and the NS moment of inertia 

Finally, one could also investigate whether there are 
other universal relations between other NS quantities. 
Recent work has shown that there is indeed a relation 
between the f- and w-modes of NS oscillations [119-121]. 
One cannot help from asking whether these relations may 
also be related to the moment of inertia, quadrupole mo- 
ment or Love number, thus yielding an I-Love-Q-f-w set 
of universal relations^". If so, one could also investigate 
whether these new quantities provide further insight into 
the fundamental reason for the existence of these univer- 
sal relations. 
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